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ABSTRACI 


A wave system extending from a body in a supersonic flow 
has a structure which is controlled by second-order cumulative 
effects and which provides for attentuation of the energy in the 
wave system and increase in the entropy of the fluid. A simili 
tude is obtained which takes into account this first-order wave 
structure; this theory appears as an extension of von Karman’s 
subsonic 


transonic similitude for supersonic flow and has no 


counterpart. This pseudotransonic similitude was first noted by 
W. Vincenti for conical flows with normal Mach Number close to 
one. As with transonic similitude, the mean-surface assumption 
is required and the theory is essentially a first-order one even 
though second-order effects are involved. The appropriate forms 
of the similitude are also given for the flow about slender bodies, 
flow at a distance from finite 


two-dimensional flow, and the 


systems 


INTRODUCTION 


oer METHODS OF analytically investigating 
the flow about bodies at supersonic speeds have 
proved themselves useful for predicting the pressure 
distributions on many types of bodies. The wave sys- 
tems predicted by linearized theory, however, have 
only geometric attenuation and do not provide for at- 
tenuation of the wave energy. Although the linearized 
wave systems may be correctly used as a basis for cal- 
culating drag they do not give even an approximately 
correct picture of the actual wave systems at an ap- 
preciable distance from the body, even though the re- 
quirements of linearized theory are clearly satisfied 
for the flow about the body. Second-order terms in the 
differential equations are needed to describe this first- 
order wave structure. 
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The classic theory of von Karman' for transonic 
similitude involves the same phenomena, but has the 
distinguishing feature that the wave structure always 
affects the pressure distribution on the body itself; 
in supersonic flow the first-order wave structure only 
affects the body significantly in the case of an edge or 
discontinuity line with the normal Mach Number close 
to one. 
may be considered to be an extension of transonic simili 


The pseudotransonic similitude here developed 


tude to cover first-order wave structure in supersonic 


flow. 


Although many results were noted earlier, particu 
larly the behavior of asymptotic N-waves, the first 
general investigation of first-order wave structure may 
be considered to be that of Friedrichs,* who gives ref- 
erences to earlier work. The theory of Friedrichs was 
extended by Lighthill,* who made a careful analysis of 
the errors of first-order wave theory in the one-dimen 
sional unsteady case; Lighthill's approach was ap 
plied to the two-dimensional steady case by Pythian.’ 
Results for the case of bodies of revolution were first 
obtained by DuMond et al.,°? and a thorough analysis 
of this case is due to Whitham.® The spherical non- 
steady case has also been studied by Whitham,’ and 
Lighthill’ has investigated the flow 
Various workers in second-order flow theory have, of 


conical case. 
course, noted certain features of first-order wave struc- 


ture. 


The present investigation stems from that of Vin- 
centi,* who noted the similitude of this paper in the 
case of a conical wing with Mach Number normal to 
the leading edge approximately equal to one. The 
approach here is one which utilizes the second-order 
equation for the velocity potential rather than the more 
generally used approach which studies the geometry of 
the characteristics and shocks involved. 
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FIRST - ORDER 


Linear and first-order wave structure 


LINE ARIZED 
Fic. 1. 


It was pointed out to the author by Vincenti after 
this paper was written that Berndt" used the pseudo 
transonic similitude parameter for the case of a swept 
wing in supersonic flow. The results of the section on 
Cylindrical Flow of this paper are thus not new but con 
firm those of Berndt. 


SYMBOLS 


al = sound velocity in the free stream 

b = wing span 

r = chord 

f, zg, h = functions describing body shape 

k = generalized similitude parameter 

s = entropy distribution function 

t = conical coordinate 

t,4 = physical cartesian coordinates 

e, = reduced cartesian coordinates 

r, 6, = reduced cylindrical coordinates 

Cy, Co, Cr = pressure, drag, and lift coefficients 

D drag 

K = similitude parameter 

V = Mach Number in the free stream 

N = cylindrical flow constant 

2 = absolute temperature and specific entropy 
U = free-stream velocity 

1R = aspect ratio 

a = angle of attack 

¥ = ratio of specific heats in an ideal gas 

n, é = coordinates for reduced wave equation 
k = additional second-order term parameter 
d = constant for conical flow 

p = fluid density in the free stream 

T = thickness ratio 

o = physical perturbation velocity potential 


yy, Vv, & = reduced potential functions 
r = thermodynamic fluid parameter 


HW, A, A = general functions for pressure, drag, and lift 


Subscripts 


normal 


nN = 
a = mean 
sl = slender body 


PRELIMINARY CONSIDERATIONS 


The disturbance due to a body in a supersonic flow 
may be considered to be composed of a superposition of 
three-dimensional wave systems. In general it is pos- 
sible to identify a reasonably small number of distinct 
wave systems, and this will be assumed to be so for the 


purposes of the considerations given here. 
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Three general principles are now stated which wij 
bear on pseudotransonic similitude. Exceptions 
these principles may be found, but their validity may }y 
accepted with minor reservations. These principles 
are : 

(1) Wave Systems Approach Planar Wave Systems 
By this statement it is meant that relations appropriat; 
to a planar wave system will be approached, not that 
the overall geometry of a wave system will approach 
that of a planar system. Focusing wave systems are to 
some extent exceptions to this principle, but the region 
of focusing action is of limited extent and the state 
ment is generally true downstream of such a region. 
Some Second-Order Effects Have a Cumulative 
Action That Is By 
meant that the action will be of magnitude comparable 
with the 
second-order effects have to be included in a first-order 


(2?) 


Unavoidable. unavoidable it is 


first-order solution. Only such cumulative 


theory. 

(3) Relative and Absolute Orientation Must Be Dis 
tinguished.—-If the geometry of a structure, in this cas 
that of a wave system or train, is affected by second 
order effects, the relative orientation or local geom 
etry of the structure is affected differently from th 
absolute orientation or location of the structure in the 
field. 

A simple example is that of two consecutive weak 
shock waves originating a fixed distance apart in a two- 
dimensional flow (Fig. 1). These shocks will eventually 
coalesce and form one shock, contrary to the behavior 
of the linear solution. The coalescence of the shocks is 
governed by cumulative second-order effects; no mat 
ter how weak the shocks are they will always come 
together. The distance to the meeting point of the 
shocks is inversely proportional to the sum of the shock 
strengths. If this pattern is crossed by a wave in the 
other characteristic direction the relative orientation is 
not disturbed. 

Another example taken from a different aerodynamic 
subject is appropriate, that of the rolling up of the 
vortex sheet behind a lifting system (Fig. 2), neglecting 
instability. The first-order structure of the 
sheet is governed by its self-induced velocities and the 
sheet rolls up as a result of the unavoidable cumula 
tive action of effects. The structural 


pattern (and relative orientation) may be shown to be a 


vortex 


second-order 


function of the reduced downstream distance 2C, bAR 
for a given spanwise lift distribution; thus, no matter 
how small the loading on the wing, there will be a dis 
tance 2 downstream for which a given pattern will ap 
pear. Absolute orientation of the vortex system de 
pends on the quantity C,/AR and would not be the 
same for wings of different loading. 

The linearized potential equation is valid for an 
arbitrary fluid, while second- or third-order forms of 
the potential equation are usually expressed for a per 
fect gas of constant specific heats. This unnecessary 
restriction will not be made here, and the quantity 
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PSEUDO 


2)(y + 1) appearing in usual second-order theories 
yd transonic similitudes is replaced by the equivalent 
yantity I’ defined for an arbitrary fluid. 

The Prandtl-Glauert 
vhich underlies essentially all similitude theories in 


coordinate transformation, 
ompressible aerodynamics, is made at the beginning 
id the analysis is carried through entirely in the trans 
jormed coordinate system. The variable = or =z is 
hosen for the axial direction. 

The rotation engendered by shock waves is of the 
third order and is properly excluded in this theory, so 
that a velocity potential exists. In this theory, as in 
the ordinary transonic similitude, the entropy rise from 
shock waves does not affect the basic equations, but 1s 
significant in that it represents the eventual energy 
lissipation and forms the basis of one method of cal 
ulating the drag. Shocks are represented by surfaces 
cross Which the pressure and the normal component of 
across which the 


the velocity are discontinuous, but 


velocity potential is continuous. 


Basic EQUATIONS 


The equation tor the perturbation velocity potential 
jor the irrotational fluid the 
stream Mach Number .1/ may be expressed, to second- 


flow of a with free 


rder terms, 


2M? 
) oss — (IV lds: = t [V'.l/°o-6-s + 
l 
ide +o (M? l)p:"f= + 
(I — 1l)dzfdes + O55 — (WP l)g:-{] (1 
where 
3 (1 /a)(Opa Op), (2) 


evaluated under free-stream conditions. The quantity 


is assumed to be always positive. In the special 
case of a perfect gas with constant specific heats the 


quantity I’ may be expressed 


[= 59 


2)(¥ + 1) (33) 


fhe potential must satisfy the condition 


¢ = 0 at infinity 


imposed upstream of the zone of action of the system of 
The boundary condition 


rv \/? — 1.2% 
go; = urs(? ’ 


( 
where 7 is the thickness ratio and c is a reference chord 


(4b 


ength, is applied for the case of a flat body, with the 
mean-surface assumption. This particular choice of 
system is not fundamental to the theory, but the mean 
surface assumption is. 

1) in Eq. (1) is really of the third 
The Prandtl-Glauert 


The term in (T° 


rder and may be dropped. 


TRANSONIC 


MILITUDE 


Sl 


coordinate transformation is made, together with a 
scale reduction by the factor c, 

\ (tV M 1) « (Sa 

y (tV M 1) « 5b 

‘ 5 /e a 


The dependent variable @ is transformed at the same 
time according to the equation 
@ = Ucr (M v 6 


and the quantities A and « are defined 


K = TM‘; (M (7 
kK (.M/ 1)/TM (S 
The resulting equation for y is 
Vir + Vw — ¥: A {pr + «ly . ve 9 
with the boundary condition 
y, = f(x, 29 on y 0 10 


PSEUDOTRANSONIC SIMILITUDE 


out ol Eq. (9) 


The limitation of the theory to those cases for which the 


No similitude comes immediately 


term in « is always much smaller than the other second 
order term would be very restrictive. The similitude 
obtained by requiring both A and « to be fixed would 
be of little practical value and could only relate certain 
systems with different values of I. 

Only a first-order theory is sought, however, and only 
such second-order effects which have a cumulative ac 
tion are to be retained. In any planar wave solution 
(v, _ y 


y.-") is zero, and for any wave system approximating a 


of the linearized equation the quantity 


planar one this quantity may be considered negigibly 
small. However, it may be of the same order of mag 
nitude as y.* in the intersection region of two wave sys 
The 


intersection with a planar wave of a planar wave system 


tems or in a markedly nonplanar wave system 


having the other characteristic direction in two-dimen 
» 


(Fig. >) 
term in y.* due to the wave itself is nonzero 


sional flow is considered as a representative 
The 
the 
clearly has 
Wy — ¥:) 


partially cancels itself if the intersecting wave system 


case. 
same magnitude all along the wave and 


is nonzero only in the intersection and 


and of 


a cumulative action. The term in 


contains both expansions and compressions; this term 
is definitely not cumulative. In general, the noncumu 
lative nature of the « term is connected with the first 
principle, that of the tendency of wave systems to ap 
proach planar ones. 

The relative orientation or structure of a wave sys 
y.*. The 
absolute orientation, however, is also affected by the 


tem is controlled by the cumulative term in 


x term and will not follow the similitude to be expressed. 
In certain exceptional cases, however, the first-order 
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wave structure may be affected by the absolute orienta 
tion; such cases are exceptions to the validity of the 
similitude. 

The «x term in Eq. (9) is now dropped* on the ground 


that it is noncumulative. The equation becomes 
Yrr + Wy — 22 = 2K 


under the boundary condition of Eq. (10). 
tions are similar with respect to first-order wave struc 
ture if, the Prandtl-Glauert 
similitude, they have the same value of A; thus the 
quantity A is the basic similitude parameter in this 
between this 


(11) 


Two solu- 


beside following usual 


theory. In view of the close relation 
theory and that of transonic flow this similitude is 
termed pseudotransonic. 
The pressure perturbation including second-order 
terms is expressed 
l T 
C, 


=G= > ve + RK? + v2 — ve) 


Me — | 2 
(12) 


k = (M? — 1)?/TM! (13) 


With the mean-surface assumption included here the 
second-order term in Eq. (12) is properly dropped in a 
first-order theory, leaving 


— (1/2)C, = r(M? — 1)7"y (14) 


For similarity under the usual Prandtl-Glauert trans- 
formation two bodies must have the same values of the 
incidence-to-thickness ratio a/r and the reduced as- 
pect ratio V 1/° — 1 AR. In pseudotransonic simili- 
tude the same must be true, leading to the general 
functional equations 


C, = r(M? — 1)~ IV M2 — 1AR, a/7, K, x, ¥, 8) 
(15a) 

Cp = 72°(M?2 — 1)~"A(V M2 — 1AR,a/1t,K) (15b) 

C, = a(M? — 1)7'?A(V M? — 1 AR, a/r, K) (150) 


Except in such cases as that of a narrow wing yawed 
at approximately the wave angle, it is to be expected 
that the parameter A has little influence upon the values 
of Cp and C,,. 

The theories of transonic and pseudotransonic simi- 
litudes have in common the facts that they are con- 
nected with the question of first-order wave structure, 
that they require the mean-surface assumption, that 
they use the Prandtl-Glauert coordinate transforma- 
tion, and that they describe the entropy increase in the 
Trefftz plane. The transonic theory is different from 


* A feature of the variational approach to linearized super- 
sonic flow is that the volume integral of (¥.2 + y,? — y,?) is 
minimized. At the New York meeting N. Rott commented that 
this approach shows the « term should have little cumulative ac- 
tion and that this might be used as part of the reasoning. 


iL 
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the pseudotransonic theory im that im the transoni 
case all wave systems are essentially merged, that the 
pressures on the body are always affected by thy 
second-order terms, and that there is no distinction be 
tween relative and absolute orientation. 

The quantity [ may be interpreted from the point of 
view of unsteady flow as the change in the signal veloc 
ity (u + a) divided by the change in uw. From this 
point of view the parameter A may be interpreted as 
the ratio of the age of the wave at a given corresponding 


point in the flow field [proportional to (¢ WV 1/* — 1) x 
(WV AT? — 1) (1/a)] and the time of transit of a 
signal velocity perturbation across the thickness of a 
wave system [proportional to (V J/* — 1/7U AM) x 


(1/T)(c/\/)]. There is an essential connection be 
tween first-order wave structure in steady and unsteady 
problems, and equivalent similitudes may readily be 
devised for various unsteady flow problems. 

In order to generalize the similitude by relaxing the 
limitation imposed by the mean surface assumption it 
is necessary to require that V 1/° — | 7 be fixed as in 
the generalized Prandtl-Glauert similitude, and thus 
that the parameter k defined in Eq. (13) be fixed. In 
this case the only nontrivial similarity obtainable is 
between flows with different values of [. In this case 
the pressure is given by Eq. (12). 

In case .\/ is very large, both « and k approach I 
and A approaches ['7.1/.. Thus, with \/ large and with 
I fixed, pseudotransonic similitude goes over into the 
correct second-order form of hypersonic similitude. 


ENTROPY DISTRIBUTION 


In linearized supersonic flow theory the drag on a 
finite system may be considered as divided into two 
parts, the vortex drag represented by the energy in the 
trailing vortex system arising from the lift distribution, 
and the wave drag represented by the energy in the 
outgoing wave system. These drags may correctly be 
computed in the linearized theory on the basis of the 
strict linearized vortex and wave patterns. These pat- 
terns are calculated at a distance away from the system 
of interest which is considered to approach infinity, 
despite the fact that with a small but finite disturbance 
the nonlinearity effectively destroys them within some 
finite distance. 

The wave system is attenuated by shock waves and 
the wave drag eventually appears as entropy rise in 
the Trefftz plane with no contribution from wave 
energy. The entropy distribution in the Trefftz plane 
is governed by pseudotransonic similitude and may be 
used to calculate the drag. 

The entropy rise across a weak shock may be ex- 
pressed 

T6S = (1/6) PW4U?[(1/2)6C, }* (16) 


which from the similitude gives the form 
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T6S r?U*Ks(x, y) (17) 
Although this expression was obtained for a single 
shock a similar one also holds for the total effect of a 


number of them, because the entropy increases may be 
17) will be considered to be the ex 


added; thus Eq. 
pression for the entropy increase in the Trefftz plane. 
The drag may be calculated by the integral 


D SS TéS di d¥ IS 


taken over the Trefftz plane. If the surface area on 
which the drag coefficient is based is equal to "AR, 


the function A of Eq. (15b) may be expressed 


A= 2QK/V IP -—1AR)S S s(x, vy) dxdy (19 


integrated over the Trefftz plane. In any case Eqs. 


(1S) and (15b) are consistent. 
SLENDER-Bopy THEORY 


Slender-body theory applies to a class of bodies of 
very low aspect ratio and of gradual contour changes, 
without discontinuities in slope. The fundamental 
characteristic of any slender body is that in the vicinity 
of the body the velocity potential is a solution of La- 
place’s equation in the two lateral variables with the 
axial variable s as a parameter. The solution depends 
for any value of z on the body shape alone with no de- 
pendence on Mach Number, apart from an additive 
function of zs which must be determined by conditions 
far from the body, which does have Mach Number de- 
pendence, and which determines the average pressure 
near the body. The Laplace equation cross-flow solu- 
tion may be considered as the superposition of line 
source, dipole, quadrupole, etc., solutions, of which the 
source distribution solution generally predominates at 
moderate distances from the axis and thus determines 
the additional function of z. The equivalent source 
distribution is given only by the cross-sectional area 
distribution; these considerations are in accord with 
the supersonic slender-body theory of Ward'! and the 
transonic equivalence law of Oswatitsch.'** The im- 
portant point here is that in most cases the linear wave 
structure, and therefore the first-order wave structure, 
is determined by the cross-sectional area distribution 
(Fig. 4 

There are some cases, as that of a very thin wing or 
of a body of constant cross-sectional area whose 
camber line changes its angle of incidence, in which the 
source solution does not necessarily always predominate 
at moderate distances from the axis. The succeeding 
term, that of the dipole distribution governing lift 
forces, will be included here. For simplicity the basic 
body shape is assumed to be symmetric with respect to 


* This law states that the solution about a slender body in 
transonic flow may be obtained from that for any other body 
with the same cross-sectional area distribution by superposing 
the difference between the harmonic cross-flow solutions for the 


two bodies 


NSONIC 


SIMILITUDE 12: 


the planes 6 = O and 6 = w 2 in cylindrical coordinates, 


about a mean camber line. This mean camber line is 
assumed to have an angle of incidence normal to the 


plane 6 = O and equal to ah(s); if the camber line is 
straight and at an angle of attack a the shape function 
h is identical to one. 

The thickness ratio is taken to be 7, The 


aspect ratio AR appears as a transverse thickness ratio 


as before 


and may in most cases be replaced by 7, because geo 
metrical conditions on the shape generally require that 
r AR be constant. In order to include the exceptional 
case of a slender flat body, however, the quantity is re 
tained, with the understanding that AR = 7 in general 
The cross-sectional area is given by 7A Rc*g(z) where 
g is a dimensionless shape function. The dipole distri 
bution will be determined in terms of the distribution 
of the area corresponding to the total effective mass 
(including the cross-sectional area itself) of the cross 
section.'* This area is expressed by A R°c*g,(z), where 
g, is another dimensionless shape function and satisfies 
the inequality g; > (7/AR)gy 

The potential equation is expressed in terms of @ in 


reduced cylindrical coordinates, in the form 
2r Af 


we oe eee ee ree 


In place of the transformation of Eq. (6) the transfor- 


mation 


@ = UcrARP(r, 6, z 2] 


is used, giving the equation 
— ©.. = 2A. PP (22 


? + (1 rib, + (1/r*)%,, 


where 


A. = TAMArAR(M — 1 2: 


is a pseudotransonic similitude parameter for slender 
bodies. 

At moderate distances from the axis the radial de 
rivative of the reduced potential may be expressed ap- 


proximately 


®P, (2 P, ( >) sin 6 
%, = a (24 
- ‘ 
The function &j) may be evaluated 
H = (1/27) go' (2 (25 


and, after some intervening analysis, ®; may be evalu 
ated 

(26 
The quantity V1/* — 1 @ AR 7 "' is an additional 
similarity parameter which determines the relative 
contribution of the dipole and source parts of the solu- 
tion in the region away from the body where wave 


systems are important. If this parameter is small it 








cs 
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may be neglected for most body shapes. The pressure 
coefficient away from the body may be expressed to 
first-order terms 


—(1/2)C, = rAR® (27) 

whence the functional equation may be given 
VM? —1aAR 
AR’ T 


C, = tAR Mh (Ky, 6.2) (28) 


with the entropy distribution in the Trefftz plane of the 
form 


T6S = 7°AR?(M? — 1)Kegissi(7, 4) (29) 


If, as has been discussed, the parameter V .J/* — 1 X 
aARr~'may be neglected, the variation with respect 
to 6 in Eqs. (28) and (29) may also be neglected. 

In the case of a slender planar body, for which the 
mean the 
t AR is small and any functional dependence on it in 
In the case 


surface assumption applies, parameter 
the formulas given above may be neglected. 
of a slender body of revolution the quantity AR may 
be replaced by 7 in the formulas above, and g; = 
2g. 

Two-DIMENSIONAL FLOW 


The case of the flow on one side of a two-dimensional 
body is now considered, with only a single wave family 
present. In this case the « term of Eq. (9) is of third 
order throughout the field and absolute orientation will 
follow the similitude; this restriction to a single wave 
family is not essential but is made for convenience. 

Dropping the lateral variable Eq. (11) may be ex- 


pressed 
Vyy — Vee = 2KvaA22 (30) 
The new variable 
= s= 9 (31) 
is introduced in place of z, giving the equation 
Vay — Wye = Bex (32) 


With only the single wave family oriented along the 
planes ¢ = constant, the term y,, is of higher order and 


may be dropped. Another change of variable is made 


n = Ky (33) 

and leads to the equation 
Woe + Ver oe (34) 
which will be termed the reduced wave equation. The 


boundary condition at upstream infinity is a null one, 
as before, and on the body is 


y, = —f(¢) on n = O (30) 


The quantity A is not a fundamental similitude 
parameter in this case, but appears as a lateral scale 


factor. The extended pseudotransonic similitude of 
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this case may be expressed: at a given Mach Number 
the product of the intensity of a two-dimensional wave 
system and the distance at which a certain first-order 
Note 
the close analogy here with the vortex system behavior 


wave structure pattern appears is a constant. 
discussed above. 


FINITE SYSTEMS 


A finite planar system is now considered, for which 
the assumption is made that nonlinear effects are un- 
important within a distance from the axis of several 
M*? — 1 


system, so that the linearized wave pattern at a dis- 


times the span or 1/\ times the chord of the 
tance from the body can be assumed to have been ap- 
proximately established. The details of this linearized 
wave pattern have been previously discussed by the 
writer in his doctoral thesis'? in connection with the 
computation of the wave drag of a finite system. Since 
the wave system has a geometric attenuation as the in- 
verse half power of the radius, it is convenient to make 
the transformation 


(36) 


y = 7 W 


in order to be working with a function which varies 
relatively slightly with r. The equation for V in cylin- 


drical coordinates is 
V.. > (1 4r7)V + (1 r-) Wee ‘on WV. = 
2Kr- "WW., (37) 
The transformation of Eq. (31) is made with r in place 
of y 
res (38) 


and leads to the equation 


V,, + (1/477) + (1/r?)Wog — 2V,- = 
2Kr~ WY, 


39) 


under a boundary condition of the form 


Vv. = -f(¢, 8) on r= 7 


where 7”) is large enough so that the linearized distant 
wave pattern has become established but not so large 
that nonlinear effects have appreciably changed the 
For r = 


(39) may be dropped as being negligible. 


r) the first three terms of Eq. 
With these 
terms dropped the boundary condition may be dis- 


wave shape. 


placed to r = O, this process probably compensating 
somewhat for the errors in ignoring nonlinear effects 


forr <r. The variable n is defined for this case 


n = 2Kr (40) 
and the resultant equation is the reduced wave equa- 
tion 

V,: + VV. = 0 (41) 
with the boundary condition 
Vv. = —f(¢, 8) on 7 = 0 (42) 
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in the form as Eqs. (34) and (35) except for the appear- 
ance of 6 as a parameter. 

The same comments are appropriate here as for the 
two-dimensional case, except that the product which is 
constant is that of the square of the intensity and the 
distance at which a given pattern appears. The func 
tion f is proportional to the radial velocity in the dis 
tant field and must be obtained using the coincident 
signal concept of the writer.'' This method may also 
be applied to bodies of revolution, without the param- 
eter 9, in which case the function f is proportional to 
the function F(y) of Whitham.’ No similitude for 
bodies of revolution other than slender ones is obtained, 
however, because in general f or F depends on the Mach 
Number and thickness ratio as well as on the body 
shape. Whitham’s method for determining the first- 
order wave structure may also be applied to the finite 
system case, but in this paper a direct characteristic 


method is used. 


REDUCED WAVE EQUATION 


The reduced wave equation, Eq. (34) or Eq. (41), 
has an immediate first integral 


Y, + (1/2)¥." = 0 (43) 


More important, however, are the characteristics de- 
termined by 


dé — ¥-dn = 0 (44) 


along which y, and y, are constant. These charac- 
teristics are straight lines in the (y, ¢) plane. 

Shock waves must be treated separately but provide 
no difficulties. The condition that yy must be con- 
tinuous may be used to show that Eq. (44) may still 
be used to describe the orientation of the shock if y, is 
replaced by the mean of its values on either side 


Ve = (1/2) (ve + Vz) 


m 


(45) 


The shock waves will, in general, be curved. With the 
boundary conditions of Eq. (35) or Eq. (42), the solu- 
tion may be readily obtained by means of a characteris- 
tic construction on the (7,¢) plane; the position and 
strength of any shock that appears is obtained by the 
solution of an ordinary differential equation. It may 
be noted that shocks can begin at a finite distance in 
the flow but that they never die out within a finite 
distance. 

The asymptotic structure for large values of 7 is 
usually an N-wave, which is bounded front and rear 
by two shock waves and has a uniform expansion re- 
gion in between. If ¢ returns to zero behind the wave, 
which is generally the case (in two-dimensions the 
body must be closed), the two shocks are of equal 
strength. The distance between the shocks varies as 
7 and hence as the square root of the intensity. In 
general the asymptotic pattern depends upon the ex- 
treme values of ¢ in the wave system. 


SIMILITUDE (27 

In the body of revolution or finite-system case the 
N-wave width varies as 7 ‘ and the shock strength as 
r‘, as was first noted by DuMond et al. 


CYLINDRICAL FLOW 


Of particular interest, because they serve to link 
pseudotransonic and transonic similitudes, are con 
siderations of cylindrical flow with the normal Mach 
Number concept. If .\/,, is the normal Mach Number 
corresponding to a given orientation on a plan form in a 


supersonic flow, the quantity 
N = [((M,? — 1)/M,?7)(M?7/(MP — 1) (46 


is preserved under a Prandtl-Glauert transformation 
Following the concepts of cylindrical flow the apparent 
thickness ratio 7, for the normal flow component is 


given by 


Mtn = Mr (47 

Defining the normal similarity parameter A, with 
respect to normal quantities, the result is obtained 

K, = TM,‘r,(M,? — 1) =N “kK (48 


and shows that A, is preserved in any Prandtl-Glauert 
transformation for which the basic parameter A is pre 
served. 

If the normal Mach Number for a given orientation 
is close to one, A, becomes the usual transonic similarity 
parameter for the normal flow. Thus pseudotransonic 
similitude corresponds to normal transonic similitude 
in cases for which normal transonic effects are present. 
including noncylindrical wings for which the normal 
Mach Number on the leading edge is close to one. In 
the case of a noncylindrical wing with the normal flow 
past a trailing edge close to sonic, the situation is dif 
ferent; normal transonic effects depend upon the ab- 
solute orientation of the corresponding wave systems. 
This absolute orientation generally does not follow 
pseudotransonic similitude where the wave systems in 
question are not in an undisturbed part of the flow 
field. Thus normal transonic effects on trailing edges 
do not follow pseudotransonic similitude in general 


CONICAL FLow 


In linearized supersonic conical flow the wave systems 
involved are of two types, of the essentially two-di- 
mensional type occurring along the characteristic 
planes tangent to the Mach cone and of the conical type 
associated with the Mach cone itself. The first type 
can mostly be handled through the two-dimensional 
and cylindrical flow developments given above; the 
second type demands special attention. If the conical 
wave system lies in the zone of action of one of the 
planar wave systems, the absolute orientation is not 
maintained under the similitude and must be related 
to the local characteristic plane orientation. Since 
this local characteristic plane orientation may be ob- 
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tamed straightforwardly and the relative orientation 
does follow the similitude, attention can be restricted 
to the case where the conical wave lies in an undisturbed 
fow. In any case, since the conical property of the 
supersonic flow about a conical body does not depend 
upon linearization, the first-order wave structure must 
be conical over the entire field. 

With the conical coordinate 


t=r/z (49) 


introduced in place of 7, and with the variable separa- 


tion 
y = 2° V(t, 6) 50) 
Eq. (11) becomes 
1— 2), + A/D — 2t? + 2nt?)¥, + 
(1 t? Wee —n(n—1)V= — 22 ‘RUN, (51) 


in which the contribution to the nonlinear term from 
waves outside the Mach cone has been neglected and 
this term has been approximated in the vicinity of ¢ = 
|. Carrying out this approximation for the linear part 


of the equation and considering @ as a parameter the 


approximate form 


(1— AW, + [In — (1/2)]¥, + 3" 'W, = 0 (52) 
is obtained, valid only near the Mach cone. Except 
for the special cases n = (1/2), (3/2) the integral may be 
expressed implicitly 

Kv, = (0) | (1 — tf) -—— —. AW, (53) 
n— 3/2 


where \ is a function of 6, corresponding to the linear 
solution, 


MOA '1 — 2) 54) 


e, = 


The case of interest here is that of regular conical flow, 
for which m = 1. In this case Eq. (53) may be solved 
explicitly 


KV, = Nw+AVA+ 1-7 (59) 


In this equation the positive sign is to be taken for a 
compression-type solution and the negative sign for an 
expansion-type solution. The latter is cut off at/ = 1, 
while the limit of the compression solution is determined 
by the shock condition 

0) 


f{—1 = (1/2)AY, 


giving the shock at 


t—1 = (3/4)\°*; AW, = (3/2)x; of) 


These solutions are shown on Fig. 5 and agree with 
those of Lighthill.. This general picture of the first- 
order structure of the conical wave system shows that 
the principal effects are of the second order in the in- 
tensity. At singular lines on the Mach cone, such as 
where a two-dimensional wave system intersects the 
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cone, a more complicated treatment is necessary which 
takes into account variations with @. The investiga- 
tion of Vincenti’ dealt with conical wings with leading 
edges near the Mach cone, for which variations in @ are 


important. 
EXAMPLES 


A) Three Two-Dimensional Airfoils 


Three symmetric two-dimensional airfoils with the 
same maximum thickness 7 are now treated using the 
One is a double 
the other 


reduced wave equation technique. 
wedge foil of chord 1 and thickness ratio 7; 
two are of chord 2 and thickness ratio 7/2, one a circu- 
lar arc foil and the other a cusped foil with two in- 
These three foils have the same 
the minimum 
Only the upper sides 


verted circular arcs. 
maximum thickness, 
values of the velocity potential. 


and hence same 


of these foils are considered. 
The first-order wave pattern for the first foil is shown 


in Fig. 6. With the foil lying between [0, —(1/2)] 


and |0, +(1/2)] on the (n, ¢) plane, the bow shock has a 
straight section from the point [0, —(1/2)] to (1, —1) 
and a parabolic section following the curve 

5S) 


n=¢ 
The stern shock is oriented symmetrically 
From the origin extends a 
fills the entire 


to infinity. 
with respect to the 7 axis. 
Prandtl-Meyer fan, which for 7 > | 
space between the shocks. 

The other two foils also have wave patterns which 
are symmetric fore and aft, and the three patterns are 
compared in Fig. 7. The solution for the second, the 
circular-are airfoil may be obtained immediately by 


considering it located between (1, —1) and (1, 1) in the 
diagram of the first. The solution for the third is 
similar to that for the first, except that a reverse 
Prandtl-Meyer fan fills the triangular region (0, —1 
(6.0), (1, —T1). 


It is easily seen from the entropy increase distribu 
tions that the drags of the second and third foils are 
equal and less than (actually two-thirds times) that of 


the first. 


B) Cusped Nose on a Body of Revolution 


In order to exemplify the slender body theory, a 


revolution at zero incidence is considered 


body of 
a cusped nose whose shape may be ex- 


which has 
pressed 


/ TV AV is.’ oY) 


in reduced coordinates. Only the flow within the zone 
of action of this tip is to be considered, and the body 
aft of the tip is assumed to have such a shape that it 
will not cause interference within the region of interest. 
The appropriate solution is given from source methods 


as 
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dD... (2 (2 — 2’)” de’ 
@ = — - Ucr : (GO) 
} a V/g'2 — 7? 
for which the approximate form near the axis is 
@ = — (5/4) Ucr?2” In (22/r) (61) 


The precise quadrature involves elliptic integrals and is 
not necessary for the purpose of this treatment. In 
the vicinity of the wave system the linearized solution 
may be expressed 


(159/32) Ucr*(s — r)?(2r) (62) 


e$<+- 
Returning to the terminology of the finite-system analy- 


sis and replacing +/V \/* — 1 by 7? in accordance 
with the slender-body theory, the boundary condition 
of Eq. (42) may be expressed 


ra 


tf = (lor/l6Y2)¢ (63) 
The solution of the reduced wave equation in this case 
is direct and essentially the same as that for the front 


half of the third airfoil of the previous example. The 
characteristics all converge to a point at 
n = 167 2/l5m (64) 
or at a value of the physical radius equal to 
cr 128 c(M? — 1) 
j ~~ 90n~2 2 (65) 
ot — 1 225 9° Pir"; 


Beyond this point a shock appears, whose position and 
strength are controlled by the body shape behind the 
cusped nose. This ¢xample is illustrated in Fig. 8. 


(C) Ring Airfoil of High Aspect Ratio 


A ring airfoil is now considered, with the usual quasi- 
cylinder assumption, and for which the chord is small 
compared with V \/* — | times the radius of the ring 
(Fig. 9). This latter restriction allows the local flow 
about the airfoil to be considered two-dimensional and 
permits the approximation of the section on finite sys- 


tems to be used. With the radius of the ring equal to 


cryo/V M2? — 1, the transformations 
v= (%9°/r*) (66) 
¢=s—f — f (67) 
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n = 2A(r* — ro’’)ro (GS 


yield the reduced wave equation with boundary cond; 
tions corresponding to the two-dimensional reduction 
The solution of the reduced wave equation gives thy 
As be 
fore, the angle 6 may be kept as a parameter, except 


first-order wave solution exterior to the ring. 


near singular points. 

For the interior of the ring the transformations may 
be altered by changing the sign of r and ro in Eq. (67 
and the sign of 7 in Eq. (68). The waves in the interior 
of the ring have a focusing action, and the approximate 
solution loses its validity near the axis. 
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ABSTRACT 


Random movements of the transition-point location on a 10 
me were observed at supersonic speeds by means of many flash 
vhlieren photographs. Distribution functions that statistically 


efine the transition-point location were determined in two 
sipersonic wind tunnels for several values of Mach Number, 
Reynolds Number per foot, and free-stream turbulence intensity 
Concomitant temperature-recovery-factor distributions were 
neasured for all test conditions 

The axial extent of the transition-point distribution function 
letermined from schlieren data increased with tunnel turbulence 
Corresponding axial extensions of the transition region de 


The data 


evel 
joted by recovery-factor measurements were observed 
suggest that a sharp transition from laminar to turbulent flow 
takes place and that this flow pattern moves randomly along the 
,.erodynamic surfaces. With the use of the statistical distribution 
functions obtained from schlieren data and surface temperatures 
is time-averaged by the surface thermocouples, a hypothesized 
instantaneous surface-temperature distribution was calculated 


for eight sets of experimental measurements 


INTRODUCTION 


pee THE STUDIES of Osborne Reynolds, many 
hundreds of experiments have been conducted to 
determine the cause and location of the transition region 
between laminar and turbulent flows. Engineering 
rules were formulated on the basis of these experiments 
both for pipes and for aerodynamic surfaces, so that 
quantities such as skin friction, which depend on the 
location of transition, could be satisfactorily predicted. 
The effectiveness of such rules often leads to false con- 
fidence among the general engineering profession that 
the transition process is well understood. Yet in spite 
of excellent pioneering research contributions by Toll- 
mien and Schlichting,' Lees and Lin,” * G. I. Taylor,‘ 
and others, the transition process has not been clearly 
defined 

Perhaps part of the difficulty has resulted from lack of 
sufficiently precise instrumentation in the transition 
region. Certainly such visual schemes as film-drying 
techniques or fluorescent lacquers could only be ex- 
pected to yield averaged information. 

In order to obtain detailed information, exhaustive 
studies utilizing hot-wire equipment appear necessary. 
Yet, with the recently expanded interest in supersonic 
flow and the almost universal application of shadow- 
the Aerodynamics Session, Twenty-Second 
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Transition-Point Fluctuations in Supersonic 


Flow 


JOHN C. EVVARD,* MAURICE TUCKER,? ann WARREN C. BURGESS, JR.7 
Lewes Flight Propulsion Laboratory, NACA 


graph and schlieren photography, observations and con 
cepts concerning transition that were not generally 
recognized heretofore have arisen without the use ot 
hot-wire equipment. 

One such observation contradicts the popular notion 
of a steady transition-point location. Schlieren and 
shadowgraph studies,®~* with short-duration exposures, 
have demonstrated that the transition point on bodies 
mounted in supersonic wind tunnels fluctuates rapidly 
While such fluctua- 


the 


with time along the body suriace. 
tions might be attributed to the turbulence of 
stream, thev have also been observed on models in 
supersonic free-flight tests conducted by Jedlicka, Wil 
kins, and Seiff of the NACA Ames laboratory. Like 
wise Schubauer (through private communication) has 
noted, by means of hot-wire pairs, the occurrence of 
transition-point fluctuations in low-speed subsonic 
flows. 

Anticipating these experimental observations are the 
theoretical studies of Emmons, and Emmons and Bry 
They visualize the transition process as re- 


son.® 1 


sulting from the growth of turbulent spots randomly 
formed in the boundary layer. The origin of these 
spots might arise from free-stream turbulence, from 
noise, or from model vibrations, as suggested in reference 
10; or possibly their origin might be traced to dust 
particles crashing at high speeds into the boundary 
layer. 

That turbulent spots can occur has been demon- 
strated by water-table experiments’ and by schlieren 
photographs of models in supersonic flow, which occa 
boundary layer 
' At least in the 


sionally show regions of turbulent 
followed downstream by laminar flow.” 
wind-tunnel experiments of Evvard, Tucker, and Bur- 
gess,®> however, the growth of isolated turbulent spots 
was seldom observed on a 10° cone in supersonic flow. 
Rather, the transition preferred to occur abruptly, 
followed downstream by flow that was predominantly 
turbulent. The not 
steady but fluctuated longitudinally along the body 


transition-point location was 
surface. 

Regardless of the mechanism whereby transition is 
initiated, the phenomenon can no longer be regarded as 
steady. Most published curves of skin friction, tem 
perature-recovery factor, or heat transfer through the 
transition region should probably be regarded as time 
averages of a complicated transition state that includes 
fluctuating transition points, multiple laminar regions, 
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or other transient phenomena. Because such complica- 
tions involve a degree of randomness, statistical studies 
of the transition region are suggested. 

The present paper includes a statistical study of 
transition-point fluctuations on a thin-walled, stainless- 
steel 10 The distribution 


functions for the transition-point locations have been 


cone at supersonic speeds. 


determined by means of a multiplicity of short-exposure 
schlieren photographs. Concurrent temperature and 
recovery-factor distributions were recorded with surface 
thermocouples mounted along a ray of the same cone. 
Eight sets of measured surface temperatures and statisti- 
cal distribution functions were utilized to calculate the 
shape of a hypothesized instantaneous surface-tem- 
perature distribution. These investigations were con- 
ducted at the NACA Lewis laboratory during the spring 
and summer of 1953. Special thanks should be extended 
to Hugo Heermann, the numerical 
calculations to obtain the instantaneous surface-tem- 


who conducted 
perature distributions. 


SYMBOLS 


g = probability density; the function g(7) dy represents the 
fraction of time that the transition point is located 
between 7 and 7 + dn 

. 
G = statistical distribution; the function G(x) = f o(n 
a 
dn represents the fraction of time that transition is 
upstream of point x 

g = average g over long time period during which flow con 
ditions change 

MI’ = Mach Number 

q = dynamic pressure 


Re = Reynolds Number 
recovery factor 
instantaneous recovery factor 


r= 
t = time 
‘ = surface temperature 

[ = average surface temperature over long time period dur- 
ing which flow conditions change 

x = radial distance from cone apex 

e@ = upstream limit of transition-point movement, ¢(7) = 0 
forn <a 

8 = downstream limit of transition-point movement, ¢(n) 
= O0forn>B 

6* = boundary-layer displacement thickness 

n = radial location from cone apex of instantaneous transi- 
tion point 

A = instantaneous surface-temperature distribution func- 


tion 


APPARATUS AND PROCEDURE 


The transition studies of this report utilized two 
thin-walled 10° stainless-steel cones. Eight sets of test 
conditions were employed for which the recovery-factor 
distributions and the statistical distribution functions 
The 
Ames cone, described in reference 11, was available for 
The Lewis cone, which was used as the 


for the transition-point location were measured. 


test 1 only. 
model for tests 2 to 8, is similar to the Ames cone as re- 
gards included angle of 10°, wall thickness of 0.032 in., 
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material of 1S-S stainless-steel alloy, and surface finish 
(maximum roughness less than 15 microin.), but js 
1S 7/,¢ in. long instead of 15 in., and internal dimensions 
differ somewhat near the apex as shown in Fig. |. The 
two cones also differ in the number and spacing of the 
constantan thermocouple wires soldered into holes jn 
the shell. A single stainless-steel wire connected to the 
model base completes the thermocouple circuit for each 
cone. 

Results obtained from either cone are believed com. 
parable. Under equilibrium conditions, the surface of 
these models closely simulated the ideal adiabatic sur- 
face in view of the stagnant air in the cone interior and 
of the thin walls and relatively poor conductivity of 
stainless steel. 

Although the instruments used with the thermo- 
couples have accuracies within +0.25°F., repeatability of 
temperature measurements during a temperature sur- 
vey was probably +0.5~F. because of gradual changes in 
tunnel stagnation temperature. The change in stagna- 
tion temperature during the schlieren observations to 
for tests 1 to 4. The 
corresponding temperature change for tests 5 to S was 
much greater, amounting to 46° for test 8S. Most of the 
variation occurred during the first quarter of the run. 

The location of the transition points for all tests was 


be described was about 3°F. 


obtained by means of a multiplicity of schlieren photo- 
graphs of about 1-microsec. exposure time. For tests | 
to 4, these photographs were taken at approximately 
uniform intervals of 20 sec. For tests 5 to 8, the time 
interval was about 7 sec. 

Tests 1 to 4 were conducted in the Lewis 1- by 1-ft. 
supersonic wind tunnel, which is a continuous-flow, 
nonreturn, variable-pressure tunnel operating at a Mach 
Number of 3.12 with the specific humidity sufficiently 
low (about 4 XK 10~° lbs. water/Ib. dry air throughout 
most of any given run) to eliminate condensation 
effects. Two cylindrical settling-chamber configura- 
tions were used to vary the stream-turbulence level. 

In the original configuration (A), the flow entered 
parallel to the settling-chamber axis, which was trans- 
verse to the tunnel-nozzle axis. A screen having a pres- 
sure drop of Sg, where g designates local stream dy- 
namic pressure, was located in the settling chamber; 
and a screen with pressure drop of 2¢ was placed at the 
exit of the settling chamber. Measurements of the 
longitudinal turbulence intensity were made at the 
Mach Number 0.12 station upstream of the throat with 
a hot-wire probe employing tungsten wire 0.0002 in. in 
diameter and 0.080 in. long. An intensity of 9.4 per 
cent was obtained for test 1. Because configuration A 
resulted in a high turbulence level, the settling chamber 
was replaced. In the revised configuration (B), the 
entrance-air direction and settling-chamber axis were 
parallel to the nozzle axis in order to achieve a straight- 
through flow. Four damping screens having an overall 
pressure drop of 10g were placed in the settling chamber. 
A honeycomb was also installed upstream of the screens. 


For ¢ 
turbt 
er U 
er C 
spect 
lock 
ypstr 
east 
f th 
sectit ) 
per Ct 
These 
tative 
vealec 
‘urve 
figura 
sugge 
Jude 
ulen 
r fro 
peak 
reflect 
rans 
The 
ttilize 
uit. 
1 SEVE 


uy 





Tes 
numb 











e finish 
but is 
enisions 
The 

of the 

oles in 
to the 


T each 


1 com- 
face of 
ic sur- 
or and 
rity of 


lermo- 
ility of 
‘e sur- 
ages in 
tagna- 
Ons to 
The 
S was 
of the 
run. 
tS was 
yhoto- 
ests | 
lately 
time 


r 1-ft. 
-flow, 
Mach 
iently 
rhout 
sation 
gura- 
level. 
tered 
rans- 
pres- 
1 dy- 
nber; 
t the 
fF the 

the 
with 
in. in 
| per 
on A 
nber 

the 
were 
ight- 
erall 
iber. 


2eT1S. 


spectively (see Table 1). 
blocks were used for both configurations A and B, the 


ection for configuration B. 
yer cent for test 2 and 1.0 per cent for tests 3 and 4. 


tative. 


figurations A and B. 
suggests that the indicated turbulence intensities in- 
dude the effects of stream disturbances other than tur- 
bulence that could arise from tunnel piping resonances 


TRANSITION-POINT 


For configuration B, measurements of the longitudinal 
turbulence intensity were also made at the Mach Num- 


yr 0.12 station, yielding intensities of 1.0 ,0.75, and 0.5 


ger cent at the pressure levels for tests 2, 3, and 4, re- 


Inasmuch as the same nozzle 


upstream turbulence intensities could be used as a 
measure of the test-section turbulence level. 


Readings 
{ the apparent intensity were also taken in the test 
These intensities were 3.5 
These apparent intensities should be regarded as quali- 
Measurements of a preliminary nature re- 
vealed the presence of a peak in the energy-spectrum 


curve below about 40 cycles per second for both con- 


The frequency range of the peak 


r from the compressor system. The amplitude of the 


peak was much lower, however, for configuration B, 
reflecting the greater effectiveness of the damping-screen 


rangement. 

The NACA Lewis 8- by 6-ft. supersonic wind tunnel, 
itilized for tests 5 to 8, has a continuous nonreturn cir- 
uit. Air is pumped through dryer beds by means of 
: seven-stage axial-flow compressor located about 150 


~I 


FLUCTUATIONS 












> SPACING 


(0) AMES CONE 


33 THERMOCOUPi F< 1” 





(b) LEWIS CONE 

Fic. 1. Cone geometries and thermocouple locations 
chamber just upstream of the throat contains a honey- 
comb and 1.8q¢ screen. The location of the model was 
approximately 30 ft. downstream of the tunnel throat. 
The Mach Number of the tunnel can be varied between 
1.5 and 2.0 by contouring the vertical flexible stainless 
steel walls. With the 17 = 1.5 wall setting, the tunnel 
may also be operated subsonically at a Mach Number 
of 0.63. 

The tunnel Reynolds Number has an approximately 
constant value of 4.3 X 10° per foot for all Mach Num- 
bers. Because there is no heat exchanger in the circuit, 
the stagnation temperature of the air stream mey vary, 



























































it. upstream of the test section. A large calming depending upon the tunnel Mach Number s-tting, the 
TABLE 1 
Test | Tunnel] Mach |Stagnation| Stagnation | Test-section} Number of | Indicated longitudinal 
number number pressure, | temperature, Reynolds schlieren turbulence intensity, 
lb/sq in. number photographs percent | 
| abs per foot | Mach | Teast 
10.12 station! section 
1 {1xl-A | 3.12 28.55 43.9 4.8x10° 183 9.4 | oe 
2 |1x1-B | 3.12 | 17.67 | 44.0 3.0 367 1.0 | 3.5 
3 |1xX1-B} 3.12 | 23.01 | 46.0 | 4-8 272 0.75 1.0 
2 [pa-B[3.12/ 401 [ 48.8 | 6.6 595 0.5 “1.0 | 
| 4 A eee ee 
5 |jexé | 0.63 | 16.8 | 161 | 3.5 
| | 17.05 | 144 | 3.6 426 
| 17.01. | | 
- + | + + —_—4 
6 |8x6 | 1.5 17.75 203.5 | 4.2 
| | 17.75 196 | 4.25 
17.97 184 |} 4.35 | 404 
17.6 174 4.38 | 
17.62 168 | 4.45 
7 [exé [1.7 20.18 161 | 4.95 | 
| 19.15 | 228 | 4.0 | 
| 19.20 | 223 | 4.1 | 360 
| | 19.26 | 215 | 4.2 | 
| 19.40 | 205 | 4.3 | 
s lexe |1.9 | 21.69 | 176 | 4.65 | 
21.69 | 226 | 4.25 
21.63 | 225 | 4.25 | 368 
21.70 | 222 | 4-28 
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Fic. 2. Schlieren of transition on cone 
ambient-air temperature, the atmospheric humidity, 
and the condition of the dryer beds. The largest 
changes in temperature occur during a given run on hot 
damp summer days, which unfortunately coincided with 
the only period that the tests described herein could be 
scheduled. 
tinent data taken during the course of the various 


Approximate operating conditions and per- 
tests are summarized in Table 1. 


RESULTS AND DISCUSSION 

The experimental distributions of recovery factor!!!” 
and local skin friction coefficients'* seemingly imply 
that the change from laminar to turbulent flow in the 
boundary layer is gradual. Inasmuch as the transition 
point has a fluctuating location®~* this implication may 
be erroneous. High-speed schlieren photographs illus- 
trate that transition is in fact a rather abrupt phe- 
nomenon (Fig. 2). 

Many photographs of the cone flow, such as Fig. 2, 
were taken at each test condition. All these schlierens 
showed abrupt transition points, the locations of which 
were distributed in a random manner in the ‘‘transition 
region.’’ From the data, the number of times transition 
occurred upstream of a given location x on the cone was 
obtained. Division of this number by the total number 
of cases (Table I) yielded the proportion of the time 
that transition occurred upstream of a specified point. 
This function G(x) is known as the statistical distribu- 
tion function pertaining to that location. A second 
statistical quantity, the probability density g(7), may 
be defined. The function g(y) dy represents the fraction 
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of the time or the chance that the transition point js 
located between y and » + dn. The relation betweey 
the two functions is given by 


G(x) = g(n) dn 


The distribution functions and the probability densj- 
ties were obtained from analysis of the schlieren data 
Temperature distributions on the surface of the cone 
were measured concurrently with the schlieren ob 
servations and will be presented in terms of recovery 
The rest of this paper is primarily concerned 
A word Ol 


factor. 
with an analysis of the data so obtained. 
caution is perhaps warranted. Transition-point loca- 
tions might be sensitive to minor changes in the cone 
surface or its geometry. Though considerable care was 
taken to prevent such changes, it is by no means certain 


the individual tests are com 


that the results from 
parable. 
EFFECT OF FREE-STREAM TURBULENCE LEVEL 


Fig. 3 compares the recovery-factor and the statistical 
distribution functions obtained at the two free-stream 
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Fic. 3. Effect of stream turbulence on transition zone 





DISTRIBUTION F UNCTION, G 
La) 


turbt 
the t 
tion 7 
facto! 
dicat: 
in the 
as rec 
fluctu 
fluctu 
turbu 
Als 
streal 
about 
Iree-S 
of tu 
porar 
ward. 
woulc 
The f. 
ous U 
transi 
streal 
The 
the sa 
on the 
the n 
Decay 
recove 
comp] 
Strate 
partly 
tions. 
higher 
streng 


E 


Dat 
show | 
distrik 
As shc 
norma 


point js 


et ween 


’ densi- 
n data 
le cone 
Cll ob 
covery 
cerned 
‘ord ol 
t loca- 
e cone 
re was 
ertain 


com- 


istical 


tream 





TRANSITION 





1.0 { _ - 
o 8 [ 
fe # 
26 
2 § 
- TEST REYNOLDS 
= NUMBER 
ry PER FOOT 
c : 6 
F 2 3.0x10 
ne 3 48 
4 6.6 
L —_—_ * SS a | 
% 6 


8 10 12 14 
RADIAL DISTANCE FROM CONE APEX, IN. 


Effect of Reynolds Number on transition distribution 


Fic. 4 
function 


turbulence levels of tests | and 3. In general, increasing 
the turbulence level increases the extent of the transi- 
tion zone about the same amount for both the recovery- 
factor data and for the distribution functions. The in- 
dications are that the rather gradual property changes 
in the transition region of such time-averaged quantities 
as recovery factor result at least in part from the random 
fluctuations of the transition point. The range of 
fluctuation is considerably smaller if the free-stream 
turbulence level is low. 

Also of interest is the observation 
stream end of the transition region (G = | 
about the same axial distance from the cone tip for both 
Certainly, random globs 


that the down- 


occurs at 


free-stream turbulence levels. 
of turbulent air passing over the model would tem- 
porarily shift the instantaneous transition point for- 
ward. Likewise, random globs of turbulent-free air 
would allow the transition point to shift downstream. 
The farthest downstream position at which instantane- 
ous transition occurs thus presumably represents the 
transition-point location for minimum instantaneous 
stream turbulence. 

The peak temperature-recovery factor occurs at about 
the same axial location as the point of maximum slope 
on the distribution function, which would correspond to 
the most probable location of the transition point. 
Decay from the peak toward the equilibrium turbulent 
recovery factor continues downstream of the G = | or 
completely turbulent point. This observation demon- 
strates that the spread of the transition region is also 
partly due to causes other than transition-point fluctua- 
The fact that the recovery factor has a peak 
turbulent value also 


tions. 
higher than the equilibrium 


strengthens this conclusion. 


EFFECT OF FREE-STREAM REYNOLDS NUMBER 


Data from tests 2, 3, and 4+ are included on Fig. 4 to 
show the influence of stream Reynolds Number on the 
distribution functions for the transition-point location. 
As shown in reference 8, the data closely approximate 
normal distribution functions for tests 3 and 4 but 
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deviate somewhat for test 2. Such deviations usually 
increase as the axial spread of the distribution functions 
increases. There is perhaps little reason to suspect that 
the distributions should be normal with such existing 
flow complications as three-dimensional model effects, 
slight tunnel pressure gradients, sound waves in the air 
stream, boundary-layer damping, and so forth. As in 
Fig. 3, the axial spread of the distribution functions 
correlated closely to that of the recovery factor. 
Increasing the stream Reynolds Number generally 
tended to the for the 
transition-point location upstream and to decrease the 
The change in width, how 


move distribution functions 
width of the transition zone. 
ever, may have been due to the change in stream-tur 
bulence level (Table |) between test 2 and tests 3 and 4. 
Based on the axial distance of the G | point, the 
Reynolds Number of transition varied almost linearly 
from 4.12 to 5.50 X 10° as the stream Reynolds Number 


was increased from 3 to 6.6 & 10°. 
EFFECT OF \MIACH NUMBER 


Inasmuch as the test-section Reynolds Number of the 
S- by 6-ft. tunnel is nearly constant for its Mach Num 
ber range, it was hoped that the data from tests 6, 7, and 
S would indicate trends of transition Reynolds Number 
with Mach Number. Such interpretations should be 
viewed with skepticism, however, because of the varia 
tions of stagnation temperature during the course of 
each test and the probable changes in free-stream turbu 
lence intensity as the compressor changed speed, pres 
sure ratio, and mass flow from Mach Number to Mach 
Number. The distribution functions for the transition 
point location are plotted, nevertheless, in Fig. 5. In 
this figure, data points close to the cone base were de 
leted, because they appeared to be influenced by the low 
base pressure. 

Presumably, if the Reynolds Number varies during 
a given test, the transition point would occur, on the 
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distributions 


average, farthest downstream at the lowest Reynolds 
Number. Utilizing these lowest values and the axial 
locations for the G = | point gives transition Reynolds 
Numbers of 6.15, 5.58, and 6.55 & 106 for Mach Num- 
bers 1.5, 1.7, and 1.9, respectively. One might presume 
from these numbers that there is little variation in 
transition Reynolds Number with Mach Number in 
this speed range. Such a presumption is somewhat con- 
trary to the trends included in references 11 and 14. 
At least part of the difference is due to the manner of 
defining the “transition point,’’ which by now is recog- 
nized as a rapidly fluctuating variable. 

At a free-stream Reynolds Number of 4.5 X& 10° per 
foot, which approximates the value of the S- by 6-ft. 
tunnel tests, the G = | transition Reynolds Number 
of the Mach Number 3.12 data is estimated at 4.6 10°. 
The value indicates an apparent decrease in the Reyn- 
olds Number of transition as the flow Mach Number 
(and many other things) changes from 1.9 to 3.12. 
This trend of decreasing transition Reynolds Number 
with increasing Mach Number above 2 was suggested to 
the authors by Lester Lees in the spring of 1952. It 
likewise conforms approximately to the trends of refer- 
ences 6, 7, 14, and 15. 

No comparison was made for the Mach Number 0.63 
data of test 5, because the G = | point was apparently 
well off the cone surface. Also, there was reasonable 
assurance from linearized theory that the model surface 
Like- 
wise, because of uncertainties in such quantities as the 
local surface Mach Number, the absolute values of re- 
covery factor should not be trusted for this test. 


was exposed to a favorable pressure gradient. 


INSTANTANEOUS TEMPERATURE-RECOVERY FACTORS 


Temperatures measured on the surface of the cone 
are time-averaged by the heat capacity of the model and 
thermocouples. These time averages are influenced by 
the shape of the probability density curve g(n). The re- 
lation between the surface temperatures and the 
probability density curve could be determined if the in- 
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stantaneous temperatures adjacent to the surface were 
known. 

For simplicity (as in reference 8) an instantaneous 
temperature distribution @ is postulated that moves 
with the transition-point location and is unaffected by 
translation of the transition point. The origin of this 
function will be taken at the transition point. Thus, 
when transition is located at a distance n from the cone 
apex, a thermocouple at a distance x from the apex is 
sub ected to the instantaneous temperature 6(x% — n) for 
g(n) dn percentage of the time. Summation over all 
positions of the transition-point location gives the time- 
averaged temperature at the thermocouple: 


*8 
T(x) = | A(x — n) g(n) dn (2) 


This average should be regarded as that caused only by 
fluctuations of the The 
limits of integration, a and 8, are chosen so that all 


transition-point location. 
nonzero values of g are included between them. 

If the stream static temperature is subtracted from 
both members of equation (2) and the resulting terms are 
divided by the difference of the stream stagnation and 
static temperatures, the time-averaged recovery factor 


results: 
*s 
r(x) = | r(x — ) g(n) dn (3) 
a 
The function 7;(x — 7) represents the instantaneous re- 


covery-factor distribution. 

Inasmuch as the time-averaged temperature 7(xv) and 
the distribution function g(y) were both measured, Eq. 
(2) represents an integral relation for evaluating the in- 
stantaneous temperature distribution 6. Before pro- 
ceeding, it may be well to note, however, that there were 
some gradual changes in free-stream Reynolds Numbers 
during the course of tests 5, 6, 7, and 8. Under such 
conditions, even though 7’ and g are time-averaged over 
a short period, they should be regarded as depending 
upon time over a long period. That is, 


bd *) 
T(x, t) = / Og(n, t) dn $) 


~“a 
Averaging both sides of eq. (4) (subject to the assump- 
tion that the function 6 does not depend on transition- 


point location) vields: 
= 8 
i(¢) = 02(n) dn (5) 


Thus the instantaneous temperature distribution may 
still be calculated even if 7 varies during the period 
for which g (n) was measured. 

Instantaneous temperature distributions were calcu- 
lated for all eight tests by numerical solution of the inte- 
gral equation. Generally, when the calculated values 
of the instantaneous distributions were reinserted into 
the equation, the calculated and measured tempera- 
tures agreed within 0.5°F. or less over the whole curve, 
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TRANSITION-POIN 


indicating the approximate accuracy of the solution. 
Nevertheless, even though in principle Eq. (2) has a 
unique solution, many instantaneous temperature dis- 
tributions can usually be found that, with a given 
probability density distribution g(m), will yield a time- 
temperature 7(x) within the 
accuracy of the temperature measurements. The 
general characteristics of the instantaneous recovery- 


averaged distribution 


factor distribution should, therefore, be obtained from 
the results of many tests rather than from an isolated 
example. 

The functions g, r, and r; are compared in Fig. 6 for 
each of the eight sets of test data. On the curves, the 
peak of the instantaneous recovery factors was arbi- 
trarily aligned with that of the time-averaged values. 

Though plots of G may appear deceivingly close to 
normal distribution functions, inspection of the proba- 
bility-density curves indicates that this was not gen- 
The instantaneous recovery factors 
Un- 


like those of reference S, the solutions for the other 


erally the case. 
for tests 3 and 4 were duplicated from reference 8. 


instantaneous recovery factors presented herein did not 
require the assumption of a normal distribution func- 
tion. 

Though some deviations are to be noted, the calcu- 
lated instantaneous recovery-factor curve is steeper and 
reaches a higher peak than that for the measured re- 
covery factor. Likewise the axial spread of the in- 
stantaneous values in the transition region is smaller. 
Like the average, the instantaneous recovery factor 
appears to approach an equilibrium turbulent-bound- 
ary-layer value at some distance downstream of the 
peak. The trend for Mach Number 0.65 (test 
notable exception. This may have resulted from the 
observation that the probability-density curve is 1n- 


Dd) isa 


complete at the base of the cone. 

The control over the tunnel flow conditions was more 
precise for tests 1 to 4, so that the resulting curves 
should be more reliable. The instantaneous recovery- 
factor distributions for these tests are superimposed in 
Fig. ¢. The preliminary conclusion might be drawn 
from this figure that either decreases in stream Reynolds 
Number or increases in stream turbulence, or both, will 
increase the axial spread of the instantaneous transition 
region. 

The occurrence of a peak on the recovery-factor dis- 
tribution raises questions concerning its cause. Perhaps 
large eddies expected to form at transition yield a dif- 
ferent proportion of heat and momentum transfer and 
a different form factor than that for a standard turbu- 
lent boundary layer, and a higher recovery factor might 
thereby result. Such a hypothesis would explain the 
relatively slow asymptotic approach toward usual equi- 
librium turbulent recovery factors. The decay process 
seems to be approximately exponential in character, re- 
quiring from 30 to 100 boundary-layer thicknesses for 


approximate completion. 


~Jj 


r FLUCTUATIONS 

A local increase in apparent recovery factor could also 
be attributed to the shock wave generated at transition 
in supersonic flow or to the tendency for the flow to 
stagnate at the subsonic flow. 
Either phenomenon would increase the local static 


transition point in 
temperature so that, if the local recovery factor were 
actually constant, the calculated recovery factor in 
terms of the free-stream conditions would increase 

Calculations were made for the Mach Number 3.12 
data to see what shock strength would be required to 
give an apparent rise of 0.02 in the instantaneous re- 
covery factor; a flow deflection of 13° would be neces- 
sary. While a change in d6*/dx at transition correspond 
ing to this deflection is not beyond reason, the long de- 
cay period is. The expansion waves generated by a 
rapidly decreasing d6*/dx downstream of the shock 
should return the flow to nearly the initial conditions 
within a few boundary-layer thicknesses of the shock. 
The conclusion persists that the change in stream con- 
ditions due to the shock is not the predominant cause 
for the rise in the instantaneous recovery factor. 

A preliminary experiment was conducted at a Mach 
Number of 3.12 to see what would happen to the re- 
covery factor if a 0.05-in. diameter O-ring were placed in 
the fully turbulent region between the peak and the 
equilibrium turbulent value. The data with and with- 
out the ring are compared in Fig. 8. 
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As expected, there was a local rise in recovery factor 
just ahead of the ring, which in this case was attributed 
to the resulting shock. Unexpected, however, was the 
large drop in recovery factor just downstream of the 
ring with a rather gradual asymptotic approach toward 
the turbulent equilibrium values. It would almost seem 
that ‘‘partly laminar’ flow was generated downstream 
of the O-ring. Conversion of turbulent flow to laminar 
has been demonstrated on cone cylinders at Mach 3 by 
Joseph Sternberg. 

A brief test was also conducted at Mach Number 0.63 
with the 0.03-in.-diameter O-ring mounted at the same 
place on the cone, but this time in a fully developed tur- 
bulent boundary layer—artificially tripped by a 0.02- 
in.-diameter O-ring placed 1'/. in. from the cone nose. 
The trend of the recovery factor near the O-ring was 
the same as in the supersonic case, but the magnitude 
of the incremental drop in recovery factor downstream 
of the ring was more than twice as great. These ob- 
servations suggest the need for additional experiments 
including heat-transfer measurements, with possible 
applications in the problems associated with aerody- 


namic heating. 


OCCURRENCE OF MULTIPLE TRANSITION POINTS 


As reported in reference 8, a study of the schlieren 
photographs from test 3 indicated that about 2 per cent 
of the schlierens (6 pictures in 272) showed the obvious 
presence of more than one transition point. The other 
98 per cent of the cases showed only one transition 
followed downstream by turbulent flow to the base of 
the cone. 

A similar study of the data obtained in the S- by 6-ft. 
tunnel indicated a more frequent occurrence of multiple 
transitions. In one case, as many as 4+ clearly defined 
transitions were observed, with alternating turbulent 
and laminar regions. A summary of the multiple 
transition observations for tests 5 to 8 is included in the 


following table: 


Number of Schlierens wita 
Total Multiple 2 3 1 
Test Mach Number of Transitions Transi- Transi- Transi 
Number Number Readings Number Percent tions tions tions 
5 0.63 126 68 16 56 12 
6 1.5 404 73 18 58 14 1 
7 b.F 360 15 12'/s 37 8 
Ss 1.9 368 55 15 53 2 


In analyzing the data, the influence of multiple transi- 
tions downstream of the first was ignored. The calcu- 
lated instantaneous recovery factors for tests 5 to 8 are 
therefore somewhat in error. 

The great difference between the frequency of 
multiple transitions in the 8- by 6- and 1- by 1-ft. tunnel 
data indicates a considerable difference in the character 
of the two air streams. As suggested in reference 4, 
small-scale eddies in the free stream would favor the 
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observation of multiple transitions if such eddies are the 
cause of the Emmons turbulent spots. A large-scale 
eddy, on the other hand, would tend to brush the whole 
boundary layer in passing across the model surface and 
would probably leave in its wake only one transition 
point. The influence of dust particles in the air stream 
is also a factor. It seems physically plausible to assume 
that the presence of dust would increase the probability 
for multiple transitions. 

This report has included considerable experimental 
evidence to demonstrate that the transition-point loca- 
tion is not steady, but rather fluctuates rapidly along 
the model surface. These fluctuations in part explain 
the spread of the transition region. 
that an instantaneous temperature distribution moves 
Whether this hypothesized 


The data suggest 


with the transition point. 
temperature distribution merely results from a mathe 
matical definition or actually conforms to a physical 
process is still open to question. The authors believe 
that the calculated instantaneous temperature dis- 
tributions do approximate a physical process but grant 
that the initial evaluations may be too crude, in that the 
influence of Reynolds Number and stream turbulence 
has been ignored. 
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Generalized Aerodynamic Forces on the 


Delta Wing with Supersonic Leading 


Kdges 


J. WALSH,* G. ZARTARIAN,? ann H. M. VOSS# 
Bureau of Aeronautics and Massachusetts Institute of Technology 


SUMMARY 


It is shown on a semiphysical basis that strip theory, that 
is, treating each chord wis¢ strip aero lynamically asa portion ol 
n infinite unswept wing, is exact in determining the contribu 
tion of any flexible or rigid body mode to the generalized aero 
dynamic force on any rigid-body mode It is also seen that this 
result is an approximation to the correct result when the second 
a flexible 


expressions are 


Employing the reverse flow 
deter 
the 


mode involved is mode 


theorem, the exact then obtained for 


mining the contribution of an arbitrary flexible mode to 
x neralized force in a second arbitrary flexible mode for a super 
sonic delta wing, with each mode represented by a double poly 
nomial series. The error involved in the use of strip theory is 
found analytically and comparative calculations are carried out 
for several cases of practical interest It is seen that the strip 


theory is sufficiently accurate for problems ordinarily encoun 


tered 
List OF SYMBOLS 

x,y = Cartesian coordinates nondimensionalized with 
respect to the chord 2/ at the root 

l = time in seconds 

w = circular frequency of oscillation in rad. /sec 

i = semichord at the root in feet 

l = free-stream velocity in ft. /sec 

V/ = free-stream Mach Number 

p = density of the fluid medium in slugs per cu.ft 

Ap(x, ¥ = pressure jump across the airfoil at point x, y 

ZalX, ¥ = vertical displacement of the airfoil at point x, 5 

0 = generalized force for mode j due to the pressure 
from mode ? 

(x, y, t) = the 7-th mode 

fi(x, y) = the 7-th normalized mode shape 

gi(t) = generalized coordinate of the 7-th mode 

\ = sweep angle of the leading edge 

gB = function defined by Eq. (12 

k = reduced frequency, k = (wb)/l 

@ = modified reduced frequency, @ = (2k.1/*) /s° 

J = Bessel function of the first kind of order » 

We = downwash distribution function defined by Eq 
(14) 

3 = the gamma function 

B = the beta function 

fy\(.\/, &) = function defined by Eq. (46 
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INTRODUCTION 


te TERM, GENERALIZED FORCE, 1S an often used, 
and at times abused, expression. Most generally 
it is the name given to the work done during a virtual 
displacement of any generalized coordinate and is most 
often encountered in minimum energy considerations 

For the purposes of this paper, only the generalized 
force due to aerodynamic effects will be considered 
The generalized force is then defined as 


x 


Oi; = 46° Api(x, vfi(x, y)dxdy 
JA. 


where Ap; is the local difference due to a 


specified deflection shape f;(x, y), and f;(x, y) 1s the mode 


pressure 


for which the generalized force is to be determined 
The most apparent application of this computation 
is in flutter analyses where the Q,,; correspond to the 
aerodynamic coupling terms in the flutter equations 
That is, the flutter mode is assumed to be of the form 


where 


g(x, y, f g(bfi(x, ¥ 3 


is the vertical displacement of the 


BAs.% 2 


is the vertical displacement of the wing in 


and 
wing, f)(v, 4 
the j-th degree of freedom, and q,(f) is the generalized 
coordinate corresponding to this degree of freedom 
In a formulation of the problem using Lagrange’s 
equation of motion, it is necessary to evaluate the 
virtual work done by the pressure distributions due 
to all mode shapes for a virtual displacement, 6q,, of 
the This 


written as 


any of generalized coordinates may be 


. 
QO; = Ap(x, y, Ofj(x, yidxdy 
1 


Ap» f,dxdy + Apsfdxdy +4 
a 1. 
l 
= (QO; 7 QV» _ @) T 
ae 


« 


Hence, each term is of the form of Eq. (1) above 


39 
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Eq. (1) may, however, also be interpreted somewhat 
more generally: let Ap;(x, y) be the pressure distribu- 
tion due to an arbitrary deformation of the wing and 
f)(x, y) be an arbitrary weighting function. Then for 
f; = constant, Q;; is the total lift due to an arbitrary 


deflection mode; for f; = x, Q;; is the pitching moment; 


for f; = y, Qi; is the rolling moment, ete. These quan- 
tities are often involved in determining the effect of 
flexibility on the aerodynamic properties of a surface. 

Throughout this work, the undeformed planform is 
taken to be thin and initially planar, and the motions 
such that the disturbance velocities are everywhere 
small compared to the free-stream velocity; that is, 
linearized supersonic theory is employed. In addi- 
tion, consideration is limited to planforms such that 
the free-stream velocity component normal to all edges 
is supersonic——the so-called supersonic planfom. While 
the strip-theory formulation holds for any supersonic 
planform, only the delta wing is studied in detail both 
for the sake of simplicity and because of its present-day 
usage. 

The evaluation of Q;; as defined by Eq. (1) for the 
delta planform with supersonic leading edges is then 
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the basis of this paper. An approximate form and its 
limitations are first found on the basis of strip theory, 
The exact form of the result is then determined for the 
case in which two arbitrary modes are involved, each 
representable by a double polynomial series. It will 
be noted that this result includes the cases of lift and 
moment for the rigid-body motions obtained previously 
by Miles.! 

Finally, the results of calculation are presented for 
several cases in which the strip-theory result is only 
approximate. 


STRIP-THEORY FORMULATION 


Consider the aerodynamic effects induced by the 
vertical motions in mode 7 of a chordwise strip of width 
dy, such as strip (AB) in Fig. 1. The aerodynamic 
forces created by strip (AB) will then lie on the area of 
the wing between the Mach lined from the leading 
edge A. It is desired to determine the contribution of 
these forces to the generalized force in mode j due to a 
mode 7 of the surface. This can be accomplished by 
determining the contribution of each spanwise strip of 
width dx, such as (CD) in Fig. 1, to the generalized force 
AB), and 


Then spanwise integra 


due to the motion of the chordwise strip 
integrating over the chord. 
tion over all chordwise strips will give the total gener 
alized force in mode j due to mode 7. 

Consider the contribution of the aerodynamic forces 
on strip CD to the generalized force due to motion of 
AB, 


(! ) (Bee ») f as ) . 
= Dp (Xo, Vo, Ts ~», wWidy 
1b? dxo e rile ’ ‘ ? . 


(5) 


where y;, and ys are the values of y at the intersections 
of the strip with the Mach lines, and Ap,'(x 
the pressure difference due only to motion of strip (AB 
of differential width. 


0» Yo, Y) 1S 


Expanding f)(%o, Vv) in a power 


series in y about vo, 


Fi(Xo, V) = Ao(Xo, Vo) + A1(Xo, Yo) (Y — Yo) + ae(Xo, Vo)(V — ve)? +... 6 


this becomes, 


l dR; i ; , ial 
= = dh Ap;‘dy + a Api’(y — ywidy + a Api'(y — yo)*dy +... (7 
4b? dX Jy J 4 ) : : : 


¥1 
But Ap,’ is symmetric in (y — yo); 
a v1 « vi 


and thus 


( 1 ) (se) 
4b? dX» 


hence, the integrals involving an odd weighting function will vanish, i.e 


” 


Api'(y — yo)dy = | Api’ (vy — w)'dy =... = 0 (S 


a f Api'dy + a { Api’ (vy — w)*dy + a | Api'(y — yo)idy +... (9 
¥1 / Vi « 


v1 
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GENERALIZED AERODYNAMIC 
fhe first integral appearing on the right-hand side of 
Eq. (9) 1s just the lift per unit length chordwise on strip 
CD) due to motion of (AB). 

But the lift on a spanwise strip of differential width 
due to arbitrary motion of a chordwise strip of differ 
ential width is just equal to the lift on the differential 
area formed by their intersection, due to the two- 
dimensional motion of an infinite straight wing under 
going the same chordwise deformation. This may be 
reasoned as follows (see Fig. 2): The pressure at point 
c due to a source, steady or unsteady, at point a is ex 
actly equal to the pressure at point 6 due to the same 
source at point d. The total pressure at c due to the 
motion of the entire strip (AB) is only affected by the 
portion of the strip forward of the Mach forecone and 
from the above is just equal to the pressure at ) due to 
motion of a similar strip at the spanwise location of 
point ¢ in the Mach forecone of 6. This is true of any 
point on the spanwise strip (CD), and hence the total 
lift on the strip due to motion of (AB) is equal to the 
lift on the differential element at b due to motion of all 
chordwise strips bounded by the Mach forecone from 
h and a straight leading edge through (A), where each 
chordwise strip has the same motion as that of (AB). 
And this is just the product of the differential area and 
the pressure due to two-dimensional motion of an 1n- 
finite straight wing. 

Hence, one may replace the first integral on the right 
hand side of Eq. (9) by the lift per unit area, Ap’'(x, 


yo), calculated from strip theory 


( | ) dR;,(Xo, Vo 
1) dx; 


A2(Xo, Me Api'(v — yvo)*dv + as(xXo, Vo) X 


« 


= o(Xo, Vo) AP™ (Xo, Yo) + 


Ap:'(y — yo) tidy + (10 


Sut do(Xo, Yo) = f)(Xo, Yo) is the deflection of the Strip 
(AB) at vo. Hence do(Xo, Yo) Ap 
bution to the generalized force which would be calcu 
lated on the basis of straight-wing strip theory. 

It may be directly concluded, therefore, that strip 
theory will give correct results for the total generalized 


(Xo, Yo) is the contri 


forces whenever the mode shape j for which the force 
is being computed is at most linear in y. Thus the 
strip theory is exact for determining the contribution 
of any mode 7, rigid or flexible, to the generalized force 
in any spanwise-rigid mode j, including roll. It would 
also seem plausible that strip theory should give a rea- 
sonable approximation to the result for modes j with 
low-order variation. This assumption is 
checked for several cases in Appendix A. 

be noted that results similar to the above 
with the limitation 


spanwise 


It may 
have been established by Miles” * 
of straight-line trailing edge normal to the stream. 
The Miles approach is wholly mathematical and differs 


considerably from that presented here. 


I 


ORCES ON THE 











PLAN 


DELTA 


APPLICATION OF STRIP-THEORY TO THE 


FORM IN HARMONIC MOTION 


From reference 4, one may deduce that the pressure 
difference in a two-dimensional unsteady flow is given, 


in nondimensional coordinates, by 


2ple™ ; e ; 
— ik e "W(x u) X 
Oo /0 

A 7 u +e ™ W(0 A “ :) : 

M M 
onion Fe 1a Fla 

) i/o “uu jdu 11 

J0 Ov M 


where Mach 
Number, w is the circular frequency of oscillation, and 


U’ is the free-stream velocity, ./ is the 


3 VM - 1 (12 
k wh l Isa 
@ = (2k\*) /B (13b 


J, is the Bessel function of the first kind of zero order, 


and I1’(x) is the downwash distribution given by 


eT UT .. 0) 
Wi(xjye™ 21k + g(x, t (14 
2b_ On 
For a chordwise mode defined by 
g(x, t Phe’ > a,(l — x (15 


i) 
it then follows that the pressure at any point along 
the chord is given by 


\ 


Zpl ~¢ s 
4px, i) = — a x 
oO 0 
| é kel —x+ u likn X 
& 0 
(l—-x+u '4o n(n —1)\1—x+u Ie 
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0,2 = we f | Ap (x, y, Ofj(x, v)dxdy (19 
A 


where f)(x, y) is given by the double summation of 
Eq. (18). 








In order to normalize the variables such that, for each 
chordwise strip, the chordwise variable is zero at the 
leading edge and unity at the trailing edge, a change of 
variable is employed, 








| n = y tan A, = (20a 


min «=p or 
i oi] Pie | 
(Il — x) = (1 — &)(1 — |n)) (20b 
It is noted that in the above the leading edge is at x = 
Q, the trailing edge at x = 1. With this change, Eqs. (17 


Consider now the generalized force according to strip 


IS), and (19) become 


. . . . - ul \ 

»¢ r > > ¢ ’ ¢ S ace > ead wt ” -\ 2 7 
the ry duc to mode 1 during N virtu il displacement of g(t, n, t) = 2be'™ SY DY anal — &)"(1 — In) 
mode j, for a delta-planform. Taking each mode as a m=0 n=0 
double polynomial series (see Fig. 3 for the coordinate ! 
system and notations), R S 

MeO y \* gi(é mt) = be YI Do bs n’ (lL — 81 = In 
; 4 P y 0 0 
gidx, ¥, t) = 2e™ 2, DL aun ee" 29 
m=0 n=0 cot A sins 


(17) . . 
QO, = 1b? cot A | | Ap; (&, 9, t) X 
JAc 


R S j 
yf y 
g(x, y, t) = 2be™ z. ys S54 = (l — x) (1S) es é on 
—, FE nn) — jn} )d&dn (23) 


0 cot A 


according to the previous section the generalized force Then, making use of Eqs. (15-16), realizing that the 
will be given exactly or approximately (depending on two-dimensional derivation is based on the local semi 
r) as chord, there results 


‘| e'"(—-4k(1 — & + a)" — dikn(l — 4 a)""'4+ n(n — 1)(1 —E4+ a)" “(Sa di + 
\Jo M 


where & and w are the values of k and @ at any spanwise station. By a slight change of variables, 
¢ = &1 — /n}),u = al — |)) (25) 


and noting that 


k = kR(1 — /n\), » = @(1 — |n)) (26) 


where & and are now the root values (vy = 0), Eq. (24) may be written in a form made useful for comparison. 


SpU2b? cot Ae’ i he "2 oe © - ntr 
QQ, = - . Zz. z, >. AmnOrs(1 — |n} — €)°n” x 
B . 1/70 m=0 n=07=05=0 
q) e~"(—4k2(1 — |y| — & + uw)” — dikn(1 — |n) —¢ 4+ 4)""' + n(n —1)1 — ni —6 + 4u)"~7] X 
J0 


we ; n if @ 
1“ u)du +(1—|n *(2ik — eae )e a (= c) fain (27) 
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EXACT FORMULATION 
The recent work of several authors, notably Ursell 
and Ward,* Flax,® and Heaslet and Spreiter,’ along with 
that of their predecessors, has led to the very useful 





y=-cotA 
principle of the reverse flow theorem in aerodynamics. n=l 
It is believed that the full usefulness of these concepts 
has not yet been realized, and, to the knowledge of the 
authors, the application to the determination of gen 
eralized forces was suggested only recently.* It is this 
ypproach which will be used below. 

The suitable form of the reverse flow theorem states, 








Apa(x, y, W(x, we“dxdy = 
JAs 
I. 


where Ap, is the pressure difference in the forward flow 


Ap,(x, ¥, OW,(x, we dxdy = (28 FIG. 4 


e 


due to a downwash IV,(x, y)e, Ap, is the pressure dif 
ference in a reverse flow over the same planform due W(x, y) . > h ( L ) a0 
to a downwash IV’,(x, ye’, and the integrals extend l 0 cot A os 
over the entire planform. This result holds for steady 


or unsteady flows, without restriction on Mach Number For simplicity in computation, the streamwise coor 


bevond those of the usual linearized theory. dinate direction may be rotated and translated, 

By comparison of Eq. (28) with Eq. (1), it 1s seen 

P = - = \ => * 
that the generalized force can be evaluated as, 


such that the origin is at the leading edge of the plan 


Qi; = 46° | | Api(x, v, Ofi(x, wi9dxdy form in the reversed flow, and the stream direction is in 
ied (29 the direction of increasing x (see Fig. 4). With this 
mm ae ¢ W(x, y change the two downwash expressions become 
= 4}? Ap;(x, y, t dxdy 
is “weoON 
erate Wwy wes 2 } 
where Ap,(x, y, t) is now the pressure distribution over hi cocoa. ne ao 
the planform oscillating at the same frequency in a re [Dik nx —"] »9 
rc —thk \ ”) 
versed flow due to a downwash distribution, 
R S ; 
W(x, y =~ y \ 
Wx, ye" = Ufx, ve (30 : = - > b : \ (34 
" 7 ? [ 0 0 cot A] 


From the modes defined by Eqs. (17) and (18), em 
ploying Eq. (14), For the flow in this direction, the plantorm is purely 
supersonic and bounded by a straight leading edge 


eee mM ooON < pe 
Wir, 9 a ae ( ) ) x normal to the flow. The pressure difference due 
. to IV’, can thus be evaluated with comparative 


[2% — nl — x) ] (1 — x)" (31) ease. 


It may be deduced directly from reference + that the pressure difference at any point will be given by 


it 
r t > » » WW (: —:- cos 4 + y) 
Ppl -e- ' ra) x 3 = | 7 
ADdx, y,t) = = Iik + i; oe *** aos vi sin @ |\dédu (35 


73 Ov 


Substituting Eq. (34), the #-integration can be carried out directly —.e., 


7 R S * 
| Regn u w ’ . 
2. p b,. ( cos 6 + ’) cos (<u sin ae oud = 


J0 r=0 0 cot A o 
R S — r - 
; ae (x — u Uu @ 
” aaa , 3 > b,. | cosé+ y cos -u sin 6 )dé@ (36) 
r=0 0 cot A J 0 Oo : M 


Making use of the identities from reference 1, 
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i [o . V r2! | é 
/ cos” 6 cos usin 6 |d6 = 1p + 1 7) 
J0 MV (* ) 2 AJ 


iu - 
av 6) 


| cos*?*! @ cos “= usin 6 |dé = O 
JW M 


it follows that 


a < w “0 fal\? _.. r! p!2° ,°* 
} cos é+ y cos usin 0) dd= 7m » ( = y - f 
Jo \pB . M p=0 \ Ba J © (2p)! (r — 2p)! T 


A 


a F , ae eae 
where B (. e+ — ) is the Beta function’ and the upper limit of summation implies the largest integer equal to or 


less than 7/2. 
Finally, employing Eq. (38) in Eq. (35) for the pressure difference and combining with II’, from Eq. (33) the ex 


act form of the generalized force is given, according to Eq. (29), by 


Spl°b? cot aa“ i : j ~ = . > ~~ m+ »p XT } 
0, =- wo pee te: ges a ee » (—n + 2Ikx) X 
B J -1J0 (m=0 Or =0 0 p=0 cot’ A 


~ 


l l 
P BY .,pt r . 
I plo! 2 2 M\?/.. 0 or a 
: ZUR + ue ie — @) ss Mu Jadu (ax ¢ 3) 
. ik 4 a I, lu pdx dy 
(2p) "(ir — 2p)! T B°a Ov] Jo MV 


wherein the change of variable 
n = ytanA 20a) 
has been made. The integrand with respect to » in Eq. (39) will be even or odd depending on whether m + r is 


even or odd; if odd, the integral vanishes, if even 


+] >! n +] *] n =] el—x 
| } = ay == 2 } | ... Gdxdy = 2 | | ... dndx 1()) 
J —le 0 J0 Jo Jd 


0 « 


and the n-integration can be carried out: 


Spl? cot Ae” “4 X & S r/2 
o,- 2 - . o.oo 
3 m=0 n=0r=0 0 rp=0 cot” A 


| tp np 5 ( ae b ] 
¢! DIZ 2 ve M 
: I (41) 
(2p) '(r — 2p)! T B*@ 


where 
>! m+) 2 l 
1 — x) P 
f= 2 | 4 x"—"(—n + 2ikx)' x 
J0 im +r — 2p + I f 


j > Py a O : sw —_ s w { ] Bis ' . - Py 
21k + ue (x nu)’ J y| —u \du rc ‘dx for m+ reven) (42) 
| ox] Jo M ‘| 


Limiting attention to (m + r) even the differentiation may be performed and the order of integration interchanged, 


f=0 (for m + r odd) 


vielding 
9 “1 fl ; 
bas da _ .)m+i 2p+lin l = , ee... ee a Ca aie 2 
’ +r—pt+i J i) iii tities |2ieu pare (i) +| 


mm ~ 
p tw I w 
suve (x — u)’ J, i 


oy 71 - 
26(s) m+? 2p+lion 1 . f 1 w 
(Il — #)” ~~ eC —e = Deenaete "SS u jdu 3) 
m+-r— 2p+ 1 Jo UU 


dxdu T 
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where 6(s) 1s defined as 


‘ . oi when s = 0 1 
10 when s + 0 
It will be noted that the x-integrations are of the form 
. *! i } 
a : a ‘ o t o j 
1 — w)%w7(x — u)'dx = | (l—wxi*(x — u Zz, (— 1] ( Vv)" edn 
; Ju t Pm (y — o)!c! f 
1)? +! *! 
= > (l — v)**" (x — u)‘dy 
a 0 (¥ _ o 'g! « ‘ 
1 (1) 7! 
= > % cng Biat+otlv+)1l-—atr 15 
o 04 —_ oO oO 


/ 


The remaining integration on uw is then of the form 


1 — u)‘u’*"e “asf uu 
Jo AM 


which can be expanded and reduced to combinations of the function 


*! - 
f,(.M, w) = | c“« a rt) 16) 
JW j 


rhe functions defined by Eq. (46) are involved in the determination of the two-dimensional air forces'’ and some 
tabulations’ |! !? as well as a ready method of evaluation!’ are available. The method of obtaining the general 
ized forces by the above reduction is illustrated in Appendix B. 


STRIP THEORY AND EXacT RESULTS 


COMPARISON | 


Consider the term of Eq. (39) for which p = 0 and call it 6 
SpU% cot de (1 PI" eh XK S | 
4 - > , oth cee n™t'x"—'(—n + ikx) X 
3 ‘ vO i 0 0 0 0 
. ? O it Ww . | = 
2k + oe uy Sol uu jducdxdyn (4 
Ox] Jo VW y 

It will now be shown that @;; is identical to Q;,°°’ from Eq. (27). If so, since Eq. (47) is exact for r 0, 1, the 


conclusion is reaffirmed that strip theory is exact for an arbitrary mode 7 which gives the pressure distribution and 
for a virtual displacement mode j which is arbitrary chordwise but no more than linear spanwise (R < | 
To prove this identity, first expand 6 
rap 9 iat M N R S 
Spb? ¢ ux a2 = - 
pl b cot Ae . N N N 


B m=0 n=07=0 0 


” -_ — 


2] >] i " 
. ° 1 Ww ] 
[ qt 'x"~'(—n + Dikx [ eeu (y e A w au | ait lin 
: 1 | J 0 VM | r { 
*] ¥ oO i 
l }n” ry" *[—n(n — L 2iknx] i} ey, — 4 a it [acy 


SpU2b? cot Ae™ > T+ = > 
8 mu=On=0r =O 0 


from which is obtained 


>! n *\ 
| ) n vr “In(n — 1) — 4tknx — 4R*x*Jeo"'(x — uy u jdugdxdn 19 
1/70 (Jo ’M { 


Inverting the order of integration, 1.e., 


*] 
h+h= | 


*\ *] . *] . 
.« dads = } ... dxdu +0 
a) 0 


« J i 
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and with 


it follows that 


1 71 n | n 
rs | / avers | (1 — [a] — 2+ w)’ 
* 1/70 J i 


Reinverting the order of integration, one finally obtains, 


| | n es 
h+d | | ‘| n”''[n(n — 1) (1 — | 
e 1/0 JW 


1k2(1 — 


Also 


| Px 
I, = | ” bia E '(—n + 2ikx) / e "(x — u) 
« l / 


>) 
= | y,” “’[—n(1 — ||)” ; 
ae 


Since /; and /; are not functions of z, z is replaced by ¢ as a symbol, and since /» is not a function of u, w is re- 


placed by ¢ as a symbol. 


Spb? cot Ae” den ae en 
eo = 3 z z z > a..h, 
A m On oy 0 0 


| bd | n vi 
(J ” iglete | a = 6 | [—4k7(1 
* 1/70 '. 0 
ni —~¢ +4)" "le (Od 


By comparison of Eq. (55) with Eq. (27 


n(n — 1) (1 


it is seen tl 
4 


Then the exact expression can be written in terms of the 
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C+ dy —<dz ay 
“n(n l)}-(1l—jn —et+u4 ‘tikn — 
s+ u)"4k?] | ¢ 1 —\n z) Yh u)és S taste 52 
L \u | 
stu “ — 4ikn(1 — |n s+ iu a 
n  —z+u)"Je l1— \|n — 3) A _ u me, oe 3 
M f 


a : ua | dn 

M ( 

| u jd Uy i4 
I t aoa in n » 


Then, combining Eqs. (53) and (54) and substituting into Eq. (48), 


likn(l — !n - ¢ + 4)’ 


n |e “ho St) bata) = 
i Hi Mw ie 7 


lat the right-hand sides are indeed identical, and thus 


n ode he ae 


(1 — ||)" | 2iR =: 


_ O 4 (56 


» strip-theory result, 


O = (),,° 4 AQ D4 


x 


where 


Splb? cot wa pe | 1 : ~: 
iy = - a }>, >, = 
* e t ( 





- ; ae ; 
> n : ay ~(—n + 2tkx) X 
0 ? cot” A 


l l 
a —, , . 
r! p!2? (; r+) M\*/.. a — @ / 
oe ere He 21k + we Se — ay se u du\dxdy (5S 
(2p) '(r — 2p)! T B°w Ov] Jo M 


Thus the error term involved in the use of strip theory 


can be evaluated explicitly. Several numerical com- 
parisons are carried out in Appendix A. 
number of cases is insufficient to be conclusive, it would 
appear that the strip-theory calculation is accurate 
enough for most applications. 

Also, although the analysis has been carried out on 
the basis of unsteady flow, it is correct for the limiting 


Although the 


case of steady flow. Thus strip theory is shown to be 
exact for such consideration as the effect of flexibility 
on the total lift, pitching moment, rolling moment 
in steady flow. In addition it may be noted that the 
reasoning behind the strip-theory formulation is also 
valid for the determination of these quantities in tran- 
sient problems, and for more arbitrary supersonic plan- 


forms 
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TABLE 1 


Comparison of the Strip-Theory Results for the Generalized Force Q;;‘* 


Ojje 
Mode Shape » £3 VU @ 2p lh? 
y 2 
20( “ Yt 3/2 O35 0.070,484,.8 
cot A 
y Pg 4 
20( = ) ew 2 O.5 0.059 227.1 
cot A 
y 2 
20( = Yi 5/2 0.5 0.049 596.4 
cot A 
ae 
2h : "3a 2 2 0.023 ,728,3 
cot A 
V 2 
20( . Yi 2 1.0 117 ,798,0 
cot A 
y 
20( = Ye 2 5 0.031 ,921,2 
cot A 
, y are 
2H ) 2 5 1. 019,747,9 
cot A 
y , , 
20( - ) (1 x)? y 5 0.003 ,838,7 
cot A 
Appendix A NUMERICAL COMPARISON OI 


Strip-THEORY AND Exact RESULTS 


The derivation above has already shown that the 
strip-theory calculation will be exact as long as the 
virtual displacement mode is at most linear in y. How- 
ever, it would also seem reasonable that the strip- 
theory calculation should serve as a good approxima- 
tion for higher order spanwise variation of the virtual 
displacement mode. In an effort to check this con- 
clusion, calculations have been carried out for several 
higher order modes. 

For the example chosen, the sweep angle, A, was 
taken as 45 In each case the quantity determined is 
Q,;, that is, the contribution of the aerodynamic forces 
due to mode 7 to the generalized force due to a virtual 
displacement in mode j. The mode considered, the 
values of the remaining parameters, and the comparison 


of Q;; and Q;; are given in Table 1. 


EXAMPLE OF EVALUATION OF THE 


Appendix B 


GENERALIZED FORCES 


with the Exact Qj; 


Percentage 


Phase Angle Oise" ™ Phase Angle Error in 
Degrees 2p lb? Degrees Magnitude 
92.3 0. 066,090 9 O18 6.2 
9] 2 0.057 648.0 91.0 ey 
GDS 0 O48 831.8 Qi) 7 15 
gO. 5 0 023 088.4 g0 4 2.4 
92.4 0.114,783,5 Q? 0 2 6 
O10 0. 080,900 , 4 GO 8 ee 
97 5 0.019, 232,6 90.7 2.6 
174.4 0.008 ,670,7 174.8 i 4 


force using strip theory would be as follows: First 
determine the generalized force on spanwise strips 
making use of previous work and available tables, then 
employ a numerical spanwise integration scheme. 
However, as indicated in Eqs. (40) to (46), the integra- 
tions in the exact expression can be reduced to known 
tabulated functions. Since the strip-theory portion 
can be distinctly identified, this will often be a simpler 
and more convenient method of evaluation. 


To illustrate, consider the generalized force Q for 


g = 20( ad ) "aig Bl 
cot A 


and (18), 


a mode given by 


Then, by comparison with Eqs. (17 
V=S=0, WV R=2 


Ao = Ay = by bo = 0, ax hy | B2 


Since the strip-theory result has been shown to cor- 
respond to the p = O-term of the exact result, from Eqs. 


11) and (43), 


SpU*b? cot Ae™ 


From the foregoing it would appear that the most 0..® = I B3a 
. ‘ “ae . ve 5 ¥: 
straightforward method of determining the generalized 8 
ind 
ik (' ['' ee a \\ ik [{' ot a 
[= — (1 — x)®<2tke "TS, u)¢ dxdu + (1 — u)Pe "So u )du 
5 Jo Ju M/s » Jo M 
tik ["' j 21k { _7(@ 
= — (1 — aw)? + l—u "ce do u )du (B3b) 
» JI 6 M 


Therefore, 
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—1( fy 
5B | 


where f, is the function defined in Eq. (46). 
When higher order Bessel functions are involved, it 
is necessary to first reduce the order, for example" 
Zip — 1) 
] (ax) = # 


aa 


J (ax) — Jy-s(ax) (B5) 


p 
In this manner, the order can be reduced to zero and 
one. Finally /; can be eliminated by partial integra- 
tion, employing either of the identities!' 


d j ' 
VSi(ax); = avJy(ax) 
dx 


(BG) 


d 


(B7) 
dx 


) Jol ax); = —aJ/(av 


Thus the integrals can always be reduced to the 


fy — function. 
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Effects of Certain Steady Motions on 
Small-Disturbance Airplane Dynamics 


MALCOLM J. ABZUG* 
Douglas Aircraft Company, Ine. 


SUMMARY 


Linearized equations of airplane motion are derived for small 
listurbances from flight with certain arbitrary steady motions 
rhe conditions studied include flight with steady coordinated 
turning rates, steady pitching or rolling velocities, and steady 
longitudinal accelerations. The equations are solved to deter 
mine the effects of these maneuvers on the small-disturbance 
ontrols-fixed modes of motion of a typical fighter-type airplane 
The main results indicate that the damping of the longitudinal 
short-period oscillation is decreased by steady rolling, but is 
relatively unaffected by the other conditions. The damping 
f the lateral oscillation (dutch-roll) is decreased by steady turn 
ing, steady positive pitching, and steady longitudinal deceleration, 
but is increased by steady rolling. Spiral mode divergence is 
xperienced in turns at high bank angles. This mode is made 
more convergent by steady positive pitching, steady rolling, and 
steady longitudinal deceleration. The rolling mode, or large 
real lateral root, is relatively unaffected by any of the conditions 
studied. These results are significant with regard to the dy 
namie stability of both piloted and automatically-controlled flight 


m extreme maneuvers 


INTRODUCTION 


geen performance standards for automati- 
cally-controlled airplanes require that flight con- 
trol systems be matched carefully to the dynamics of 
the airplane. The classic equations of airplane motion 
for small disturbances from steady straight symmetric 
flight' are generally used in system design to calculate 
the modes of airplane motion and airplane transfer func- 
tions.” * However, for airplanes which must be auto- 
matically controlled during extreme maneuvers, a sta 
bilization system design based on small disturbances 
from steady straight symmetric flight may be inade- 
quate. At the same time, it is desirable to avoid the 


difficulties inherent in using the general nonlinear 


equations of airplane motion for this purpose. 

In this paper, linearized equations of small-disturb 
ance motion are developed which are applicable to small 
disturbances from four conditions of maneuvering 
flight. The effects of 
trols-fixed modes of disturbed motion are investigated 


these maneuvers on the con- 


for a typical fighter-type airplane. 

General small-disturbance equations of motion were 
previously developed by Frazer, Duncan, and Collar, 
using matrix techniques. However, those equations 
involved disturbance variables in terms of earth co- 
ordinates, and are considered of less general usefulness 


than those of the present paper. 


Received January 14, 1954 
* Design Specialist, Aerodynamics, El Segundo, Calif. 


SYMBOLS 


NACA stability axes® are used in this analysis. Vector quan 


tities are indicated by boldface type 


aerodynamic force vector 
acceleration of gravity vector 
angular momentum vector 
aerodynamic moment vector 
total velocity vector 
disturbance rotation vector 


angular velocity vector 


KO xO WR 


multiplication sign for vector cross-product 

D absolute change (seen on the earth reference 
frame or on a reference frame in steady 
linear motion 

relative change, rate of change (seen on the 
airplane reference frame 

acceleration of gravity 


moment and product of inertia about the 


\ P ) ° Z and A Z axes 

UL. MN rolling, pitching, and vawing aerodynamic 
moments 

ni mass 

\ 0.V/Op, etc 

nN load factor 

r,t & rolling, pitching, and vawing velocities 

Dp, g,? disturbance rolling, pitching, and vawing 
velocities 

operator d/dt is the Operator Sat 

u,v, U disturbance velocities along the Y-, }-, and 
Zz axes 

V total velocity 

On Fe 4 longitudinal, lateral, and normal aerody 
namic forces 

a flight path angle 

o bank angle (rotation about Y-axis 


Subscript 


] steady flight condition about which dis 


turbances occur 
DERIVATION OF EQUATIONS 


General Equations of Airplane Motion 


General equations of airplane motion which apply at 
any instant’ ° may be written in vector notation as 


F + mg = m(dV/dt) +2 XK mV\ 


Q=dH/it+QxH § !) 
General Small-Disturbance Equations 
In disturbed flight, 
F = F, + DF 
&é =e + Dg=& (2) 
V= V, + DV, etc. 
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TABLE 1 - 
General Equations of Small-Disturbance Airplane Motion. 
COEFFICIENT OF 
Sane aes ae ere ee —_ pee 
“u vu w P F ” 
=4X An ios 
, s -R, Q, fe) cos Si cos , s - v % 
+ (4Vdt) = F > g cos &, SIN ,s 
on -! 
S4v/m = R, s -P, - g.cos ¥, cos $s ) v,- [qsin 5+(av/at) js 
rT 
a e “1 
“ 24 Z/m= -Q, P, s g Cos 8 SIN > s |-v, +(g sin 5 +(avfae)}s fe) 
z 
fe) ee a ea — 
be 
5 | a 
a Sal = 0 0 ) 3,8 -°@,Ty, | R(Lz-Ty)-F Lyz -I,,8+9 (I, -I,) 
lJ 
| 
Sam = 0 ) ) R,(I,-1,)t2PAIy2 | 3ys PB (I,-T,)-2R, Ty. 
SaN = fe) fe) O -1,,8+@,(I,-I,y) P(I-Iy)+R,Iyz Iz7S+Q,Iyz 





























where F, g, V, etc., apply to disturbed flight, F,, gi, Vi, 
etc., apply to the steady flight condition that existed 
previously and DF, Dg, DV, etc., are the absolute vec- 
tor changes that occurred in the interval. The abso- 
lute vector changes are related to changes seen by an 
observer on the airplane reference frame dF, dg, dV, 
etc., by :§ 


DF dF +X X F, ] 
Dg dg +X X gi (3) 
DV =dV+xX xX Wi, - 


II 


Substituting Eqs. (2) and (3) into (1), neglecting 
the products of infinitesimals, and simplifying, the gen- 
eral equations of motion for small disturbances from 
steady flight are: 


dF/m + X X [(dV/dt),; — gi] + 
(dV/dt),; = dV/dt + Q; X dV + 
(4) 
dQ xX ™ 
dQ = dH/dt + Q, X dH + dQ X Hi, 


Eqs. (4) are given in Cartesian form in Table 1, using 
NACA stability axes.* These equations apply, in 
general, to small disturbances from asymmetric flight 
since the initial rolling and yawing angular velocities 
are not necessarily zero. Asa result, the aerodynamic 
cross-derivatives such as /, and M, are also not 
necessarily identical to zero. 


Classic Small-Disturbance Equations of Airplane Motion 


The classic equations of airplane motion for small 
disturbances from steady straight symmetric flight 
may be obtained from Table | by setting the following 


quantities equal to zero: P, Qi), Ri, gi, and (dV /dt),. 


In this case, the equations are separable into symmetric 
or longitudinal, and asymmetric or lateral sets. 


Equations for Small Disturbances from Turning Flight 
In steady turning flight about a vertical space axis: 


P, = —Qsin ¥1 
Q, = 2Qcos y; sin gd (5) 
R, = Qcos y; cos 6) 


where, for coordinated turning flight (resultant aero- 
dynamic force in the plane of symmetry), 2; may be 
approximated by: 


Q, = 7 tan ¢; V; (6) 


Eq. (6) is valid for large absolute values of the initial 
flight path angle 7. 


Equations for Small Disturbances from Steady Pitching 
or Rolling Flight 


These equations may be obtained from Table 1 by the 
appropriate substitutions of Q, or P;. In the former 
case, the equations are separable into longitudinal 
and lateral equations. The equations of lateral motion 
for small disturbances from steady pitching flight 
may be used to obtain a linearized representation of 
the rolling-pullout maneuver. The airplane’s load 
factor is found in this case as: 


n = Q:Vi/g + cos y (7) 


A simplified form of the equations of motion for small 
disturbances from steady rolling flight has been de- 
veloped by Phillips.* 
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STEADY BANK ANGLE,®,, DEG. 


Effect of steady turning on modes of motion of typical 
fighter-type airplane 


Fic. 1. 


Equations for Small Disturbances from Longitudinally 
Accelerated Flight 
These equations are obtained from Table | by letting 
P,, Q1, Ri, and ¢; equal zero. The equations are then 
separable into longitudinal and lateral equations. 


RESULTS 

The effects of the conditions studied on the controls- 
fixed modes of motion of a typical fighter-type airplane 
have been determined by making numerical solutions 
of the equations of Table 1. These calculations were 
made assuming zero airspeed change during the mo- 
tion, for three flight conditions. Stability derivative 
values used in the calculations are listed in the Ap- 


pendix. 


Effects of Steady Turning on Modes of Motion 
The effects of steady turning on the modes of airplane 
motion are presented in Fig. 1, as functions of the steady 


CERTAIN 
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bank angle. The lateral oscillatory mode becomes 
more lightly damped, and increases in f.equency as the 
bank angle is increased, while the spiral mode becomes 
slightly more convergent at first, and then divergent at 
extreme bank angles. The other modes of motion are 
substantially unaffected, except for the long period or 
phugoid oscillation. Although not shown, the damp- 
ing of the phugoid mode is sharply increased at high 
bank angles, due to the effect of increased drag. 

An important result of this study is that the lateral 
and longitudinal modes of motion are not independent 
for small disturbances from steady turning flight. For 
example, the lateral oscillation following a_ sideslip 
angle disturbance affects measurably the angle of attack 
and pitching velocity variables at a high steady angle 
of bank. This coupling is evident in Table | in a num- 
ber of terms involving both gravity forces and inertial 
forces and moments. The characteristic equation is a 
ninth-degree polynomial, with a single zero root repre- 
senting the neutral azimuth mode. 
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Effect of steady pitching on modes of motion of typical 
fighter-type airplane. 
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Effect of steady rolling on modes of motion of typical 
fighter-type airplane 


Fic. 3. 


Effects of Steady Pitching on Modes of Motion 


The effects of steady positive pitching on the modes 
of motion are presented in Fig. 2 for the example air- 
plane, as functions of the steady pitching velocity. 
At constant airspeed, the short-period longitudinal 
Slight reduction in damping of 
the lateral oscillation is shown in Fig. 2 for positive 
due to a negative contribution to V,. On 


mode is not affected. 


pitching, 
the other hand, the spiral mode becomes more conver- 
gent due to a negative contribution to the derivative 
L,, for the particular inertia values considered. The 
rolling convergence is relatively unaffected. 


Effects of Steady Rolling on Modes of Motion 


The effects of steady rolling on the oscillatory modes 
of motion shown in Fig. 3 agree substantially with the 
simplified calculations of Phillips,* for the same con- 
stants. As indicated by Phillips, the higher-frequency 
mode (longitudinal) increases in frequency at high 
rolling velocities, and conversely for the lower-frequency 
The reduction in damping of the 


mode (lateral). 
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longitudinal oscillation shown could be serious for g 
highly maneuverable airplane. Steady rolling makes 
the spiral mode strongly convergent, for the low-speed 
flight condition. As in the case of small disturbances 
from steady turns, the longitudinal and lateral modes 


of motion are coupled during steady rolls. 


Effects of Steady Longitudinal Acceleration on Modes of 
Motion 


The effects of steady longitudinal acceleration on the 
modes of motion of the example airplane are shown in 
Fig. 4. The reduction in damping of the lateral oscil- 
lation shown for deceleration may be a problem on air- 
planes where sudden reductions in thrust or increases 
in drag are possible. For the same conditions, the 
convergence of the spiral mode is greatly increased by 
an outward-acting effective side force due to yaw. 


CONCLUDING REMARKS 


The effects of steady turning, steady pitching or roll- 
ing velocities, and steady longitudinal acceleration on 
the small-disturbance controls-fixed modes of motion 
of a typical fighter-type airplane have been investi- 
gated. From the standpoint of both piloted and auto- 
matically-controlled flight the most significant effects 
found were: 


(Continued on page 762) 
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JOHN W. MILES' 


Unwersity of California, Los Angeles 


ABSTRACT 


Experience has shown that the fluctuating loads induced by a 
jet may cause fatigue failure of aircraft structural components 
In order to throw some light on this and similar problems, the 
stress spectrum and the ‘‘equivalent fatigue stress’ of an elastic 
structure subjected to random loading are studied. The analysis 
is simplified by assuming the structure to have only a single de- 
gree of freedom and by using the concept of cumulative damage, 
the results being expressed in terms of quantities that can be 
directly measured. As an example, a similarity expression for 
the probable value of the equivalent fatigue stress of a panel 
subjected to jet buffeting is derived 


NOTATION 


! = panel area 
1 (w) = resonant amplification factor, Eq. (2-12 

"i = undetermined coefficients 

D = cumulative damage, Eq. (4-3) 

E = Young's modulus for panel 

F(t = instantaneous magnitude of force 

Fy = root mean square load—i.e., V Fx) 

V/ = subscript referring to Miner (reference 12) 

Vis = number of complete stress reversals of amplitude s 
required to produce fatigue failure 

A = Vs 

Pi 4 the probability of the event | } 

‘a = probability distribution of the envelope of y(t), Eq 
3-4) 

Re = Reynolds Number of jet 

/ = period of dominant mode of oscillation 

l = velocity of jet with respect to surrounding medium 

Z = impedance of oscillator, Eq. (2-11 

= sonic velocity in air surrounding jet 

d = jet diameter at nozzle exit 

erfe = complementary error function 

{(w = power spectral density of the force F(t 

falw = see Eq. (2-8) and Fig. 4 

fp(w) = see Eq. (2-9) and Fig. 4 

¢ = gravitational acceleration 

e(w power spectral density of y(t 

h = panel thickness 

k = a constant defined by Eq. (4-6); & = 1 on Miner's 
theory (reference 12); k = 2 on Shanley’s theory 
(reference 13) 

n = number of complete stress reversals at amplitude 

p = static pressure of air into which jet exhausts 


stress 


= instantaneous stress 


Presented at the Aeroelasticity Session, Twenty-Second Annual 
Meeting, IAS, New York, January 25-29, 1954 

* The present paper differs in several respects from IAS Pre 
No. 435 and Douglas Aircraft Report SM-14795 (June, 
both of which were titled ‘‘An Approach to the Buffeting 


print 
1953 ), 
of Aircraft Structures by Jets.” 

+ Associate Professor of Engineering; also Consultant, Doug- 
las Aircraft Company, Inc., Santa Monica, Calif 


= equivalent or reduced fatigue stress; see Section 4 
= static stress produced by F 
s = stress level (hypothetical) at which fatigue failure 


occurs in one complete cycle 
time 
= specific weight of panel material 
ierodynam 


exhaust) of 


x = distance from source (jet 
noise to panel 

yt = instantaneous displacement of oscillator 

y = static displacement that would be produced by the 
root mean square load ( F 

y(t = mean square displacement 

I = Gamma function; (nm + 1) = nm! if ” is an integer 

a = see Eq. (4-1) and Fig. 8 

6 = damping of dominant mode expressed as a fraction 
of critical damping and equal to logarithmic dec 
rement divided by 27 (if 6 l 

p = mass density of air 

> = an arbitrary function 

¥ = mean square value of displacement or stress 

w = angular frequency 

w = angular frequency of dominant mode of oscillation 

w = characteristic frequency of input spectrum f(w 


equal to frequency of fmar. unless fmar. at w = 0, 


in which case w, designates (approximately) half 


9 7) 
«o~é 


power point, Eq. ( 


(1) INTRODUCTION 


tbs PROBLEM OF DETERMINING the stresses resulting 
from random loading and the consequent possi 
bility of fatigue failure appears to have received but 
However, recent 
the 


little attention in structural design. 
advances in both the statistical techniques’ * and 
basic physical concepts*~* render possible an approxi 
mate analysis of some of the more important features of 
the problem. The statistical techniques, at least in 
their modern form, are due primarily to Rice! and Wie 
ner.” The physical concepts underlying the phenomena 
to be considered have recently been reviewed and ex 
pounded by Liepmann* * and by Luskin and Lapin 
the present work may be regarded as an extension of 
that of Liepmann.* Reference also may be made to 
the fundamental work of Lin,® Uhlenbeck and his stu 
In addition to 
these references, the writer is indebted to H. Klein, 
H. W. Liepmann, and H. Luskin for frequent discus 
sions of this problem. 

The present 
possibility of panel failure under the fluctuating pres 
sure (i.e., buffeting) of a jet exhaust (Fig. 1), a special 


dents,’ and Graham and Rodriguez.” 


investigation was prompted by the 


case of the more general problem of structural fatigue 


under random loading. A second example of aero 
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4F(t)= RANDOM FORCE RANDOM FORCE 


oA Are = 
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t= TIME PANEL 


Fic. 1. A panel subjected to a random force resulting from 
the fluctuating pressure field of an adjacent jet. In the analysis 
of the text it is assumed that the force may be considered to act 
ata single point 


nautical interest is repeated gust loading of a wing. 
The general problem can be resolved into three more 
or less distinct parts: 

(a) The statistical description of the random loading. 

(b) The statistical description of the dynamical 
response of the prescribed structure under this loading. 

(c) The determination of the probability of fatigue 
failure associated with this response. 

Problem (a) will not be treated in the present study, 
although two representative statistical descriptions of 
the loading will be considered. In general, this in- 
formation must be obtained experimentally (Fig. 2), 
but theoretical considerations also may be valuable. 
If the excitation is aerodynamic noise (as for the jet 
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SPECTRUM 
‘cecil F(t)) ANALYSER 
AF*: f(a) Aw 
RANDOM 
INPUT f(w)=POWER SPECTRAL 
DENSITY 
W=FREQUENCY 
@ 
Fic. 2. The experimental determination of the power spec 


tral density f(w), of a random function, F(t). F(t) is recorded 
(e.g., on a tape) and fed into a spectrum analyzer that transmits 
only those components within the pass band @ + (1/2)Aw. The 
intensity (to which the power is proportional) of the out- 
put in the pass band—viz., AF? 
meter, and the power spectral 
by dividing AF? by Aw. Thus, the meter reading, AF?, appears in 
the power spectral density vs. frequency plane as an increment 
of area of width Aw and mean ordinate f(a). 


is recorded on a mean square 


density at @ is determined 
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exhaust already cited), there exists the recent theory 
of Lighthill,? which has received some experimental 
verification.'” Gust loading and related problems have 
been discussed extensively in references 3, 4, 5, and 1] 
and also in the references cited therein. 

The determination of the response of a dynamical 
system under random loading has been considered in 
The for the 
study is a single degree of freedom (lightly damped) 


references 3-S. model selected present 
oscillator loaded at one point, for which the analysis 


This model 


is similar to that of references 3 and 6. 
exhibits the essential feature that distinguishes the 
dynamic response of an elastic system from the corre- 
sponding static response—namely, the relatively sharp 


(frequency) selectivity associated with resonance. 
On the other hand, a distributed structure has an in- 
finite set of resonant frequencies and a (possibly sev- 
eral) characteristic length(s) that must be compared 
with some characteristic length of the excitation in 
order to establish a measure of the phase differences 
The 


coming of the single degree of freedom model is not 


among loads at various points. former short- 
believed to be serious, since one mode usually will 
dominate the dynamic response of a distributed struc- 
ture, but the effects of distributed loading would ap- 
pear to justify further investigation. 

At the present time, the description of fatigue appears 
to rest on an almost entirely phenomenological basis, 
The 
present analysis will be based on Miner’s concept of 


and analytical work on the subject is slight. 


“cumulative damage,’''? and on Shanley’s extension!’ 
of this concept. It is the writer’s opinion (cf. refer- 
ences 12, 13, 17, 19, and 20) that the cumulative dam- 
age concept is probably adequate in the sense that 
presently available data do not permit the prediction 
of fatigue life to better than a factor of two (indeed on 
statistical grounds, this may never be possible for a 
single sample); however, an alternative approach is 
discussed in the Appendix, and it is concluded that the 
assumption of a representative type of nonlinear 
accumulation of damage leads to only an insignificant 
change in fatigue life when, as in the present study, 
stress amplitudes are distributed continuously over a 
wide range. 


(2) RESPONSE OF A SIMPLE OSCILLATOR 


A simple oscillator can be characterized by wy», its 
resonant frequency, 6, the ratio of the actual damping 
to critical damping (and equal to the logarithmic dec- 
rement of the free oscillations divided by 2 7 provided 
that this decrement is small), and yo, the static response 
to some reference input Fo. If the mechanical model 
of Fig. 3, having mass m, spring constant k, and vis- 
cous damping constant c, is adopted, then 


wy = VR/m 


6 = c/2V km (2-2 
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haary ; If, on the other hand, f(w) exhibits a peak at some 
COTY PAA LLLLLLLLLLL ' ae : ‘ 
Nental frequency w,, that is appreciably greater than zero it 
> have is convenient to choose w; = w,,. 
RESONANT ; ; ; 
nd 1] C2, = FREQUENCY Examples of these two possible types of power spec- 
K C tral density are given by 
mucal oe (a) Monotonic: 
ed in 
‘esent f.(w) = (2F?/V «ewe ”* (2-S 
aped) S _ RELATIVE 
alysis DAMPING (b) Peaked: 
nodel ee 
: tl I,(w) = (4F? Vv 7Wa;) (Ww w)" e wits (2-9 
) 1€ m 
Orre- — ui , eee 
; . and are plotted in Fig. 4. The half power point for 
shar 5 i fa 
P ‘= +): STATIC f, occurs approximately at w = 0.48 w. It is, of course, 
ance, ” ine ' 
y E/ RATIO conceivable that /(w) could exhibit several maxima, 
nN in- 
but f,(w) and f,(w) may be considered to be adequate 
sev- . ar - 
' a representations (for the present investigation) of the 
are ; ; 
: oa Fic. 3. A model of a single degree of freedom oscillator, char more important types. 
n in acterized by its resonant frequency, relative damping, and static os : , ‘ ? 
tio The response of the oscillator at a particular fre- 
nces os ; ee hie ' ; 
at quency w is governed by its force-displacement im- 
ion pedance,* Z(w), and the power spectral density of this 
— ae 9.2 : R . b A 
will Yo = Fo/k ‘--?) response is given by (reference 3) 
ruc- but, irrespective of the actual model, the parameters : ‘ 
: : A Bi) ie g(w) = f(w)/| Z(w) (2-10) 
ap- w), 6 and (yo/ Fo), being more amenable to generalization 
and subject to experimental determination, are prel- |) terms of the parameters defined above, | Z|? is given 
ears erable to such a set as (m, c, R). ber 
iSiS, Let y(t) represent the displacement of the oscillator 
The from its position of static equilibrium (mg/k for the , Fy\?4 w \272 w \*) 
t of model of Fig. 3); then the differential equation govern- Z(w)|? = (" ee: (“ + 46° ss ) ( 
mn! ing its response to the input F(t) is, in dimensionless . ; (2-11 
fer- form, 2 
oe F ‘ ‘ Fit) The ‘‘power amplification curve,’ or the square of the 
lat : + 3i-— 4-5? — (2-4) usual resonance curve, is given by 
ion Wo" Vo wWoVo Vo Fy 
on = e ° wae ) : ° x 4 . “J . ‘e¢ —_ ,twi ae wet 7/4 
which, for the model of Fig. 3, is equivalent to the more Z (w) is defined such that if F(t) = e'', y(t) = e''/Z(w 
yr a ° eq: 
: familiar form 
is 
the my + civ + ky = F(t) (2-5) 
“ar ; ; : - ; fiw) = POWER SPECTRAL DENSITY OF 
- However, the dimensionless form, Eq. (2-4), is to be 7 INPUT FORCE 
™ preferred in the present study, not only in being inde- > 
ee pendent of the particular model but also in that, by Ft (a) MONOTONIC - 
virtue of Hooke’s law, the displacement y(/) may be , MEAN SQUARE: F#/tiard 
replaced by the stress s(t) provided that yo also is re- oi 
I : I ; 
placed by so, the static stress produced by F). b) RE AKED 
Now suppose that F(/) is a random function having (b) Ps 
ts the power spectral density f(w). The latter can be 
ig broadly characterized by two parameters: the mean 
.- square value of F(t), given by 
¢ | (t), g j AF* 
d — F 
- Fxt) = J,” fw) dw (2-6) 
e] a : 
and some characteristic frequency, here designated as 0 i Z > 
az a ° . . e . » Y 
oe “ally decreasing fanc O w,= CHARACTERISTIC w = FREQUENC 
oe If f(w) is 1 monotonically decreasing unction of FREQUENCY 
w, it is convenient to choose w; = w:;, Where w.,, defines 
(at least roughly) the ‘‘half power point’’—viz. Fic. 4. Typical power spectral densities, as given by Eqs 
F (2-8) and (2-9). Such curves usually are determined expert- 
=) Oe . mentally, a single reading of the type illustrated in Fig. 2 lead- 
} So J(w) dw = Se , I(w) dw (2-7) ing to an increment of area such as that shown under curve (b 
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The power amplification curve of a single degree of free 
dom oscillator, as given by Eq. (2-12). 


. 2-12 
{1 — (@ wy)? ]? + 16-(w ea 
and is plotted in Fig. 5 for 6 = 0.1. 
If (6w) w,) = O(1), as illustrated in Fig. 6, the re- 


sponse of the oscillator is essentially resonant and it 
receives appreciable loading only in the neighborhood 
of resonance, in which neighborhood its power spectral 


density is approximated by 


(ax) Vo 4 6w 
r _ * f(ay)A(w), = ((1) 
Z(w) ; Fy w] 


(2-13) 


~~ 
U ’ — 
£(w) 


and (in the neighborhood of o ) has the same shape as 
the (power) resonance curve of Fig. 5. If, on the other 
hand, (6w) w;) >> 1, the response of the oscillator would 
be essentially static (4 = 1), and its power spectral 


density would be given approximately by 

f(w) oe | 6a 

r _ : I(w), Pe 
Z(O) Fy w) 


and would have the same shape as the power spectral 


~ 
g o) sum 
( ) 


density of the loading. 

The first of these two possibilities is technically the 
more interesting, since the amplification associated with 
resonance is generally large. Thus, resonance appears 
to be an essential factor in the fatigue failure of a panel 
subjected to jet buffeting. On the other hand, the 
loading of a wing by randomly encountered gusts 
might be considered to be essentially static (but the 
motion of the airplane as a rigid body may be impor- 
tant) in consequence of the relatively low frequency 
(compared with, say, the fundamental wing bending 
frequency) at which gusts of appreciable magnitude 


are met. 
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The subsequent analysis is directed primarily toward 
the problem of essentially resonant response, since it is 
not only of greater practical interest but also is mathe- 
matically simpler by virtue of the concentration of the 
power spectral density of the response in the neigh 
borhood of the resonant frequency. 


(3) STATISTICAL PROPERTIES OF THE RESPONSE 


The complete statistical description of a random 
function y(‘) requires the specification not only of its 
power spectral density g(w), giving its frequency dis- 
tribution, but also its probability distribution in time. 
The distribution most commonly encountered and also 
most amenable to analysis (c/. references | and 2) is 
designated as normal (or Gaussian) and gives as the 
probability that y(t) may be found between the pre- 
scribed values y and y + dy 


' 2yo7 ‘ 
Piv<yvt)h<yt+tdy =e dv NN 2ry 


where Y designates the mean square value of y(/ 


ViZ., 


te = 9) = So"ele) de 3.2 


Certain results of practical interest may be caleu- 
lated directly from Eq. (3-1) without the necessity of 


further assumptions. Thus, replacing the displace- 
ment y by the stress s and regarding y 


square stress, the probability of static failure, defined 


as the mean 
by s(t); > s,, where s; is some arbitrarily prescribed 
value (e.g., elastic limit or yield point), is obtained 


by integrating Eq. (3-1)—-viz., 


>! r ' )» j € wd 9 
P} S(t) >S = 2 y-oa) 


However, the calculation of such quantities as probable 
fatigue life demands that additional information and or 
assumptions be introduced. 

The analysis of the essentially resonant response is 
especially simple in that it is completely characterized 
by the resonant frequency w) and the mean square re- 
sponse Y (cf. the broad characterization of F(t) by 
w, and F? Moreover, whatever the statistical dis- 
tribution (in time, not frequency) of the input /(f), the 
response y(/) tends to a normal distribution as its band- 
width becomes small compared with that of f(w)—1i.e., 


as (6a) w,) tends to zero.* 


* It is generally argued that this result follows directly from 
the central limit theorem (see also reference 14). 
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a. 
® 


Fic. 6. Illustrating the approximation, Eq 
the relatively narrow bandwidth of the oscillator amplification 
factor compared with the bandwidth of the input [assumed as 

w) in the sketch, but the result would be essentially similar 
for an input spectrum of the monotonic type / 


A qualitative picture of the essentially resonant re- 
sponse to a random input may be formed by noting 
that a lightly damped oscillator acts as a narrow band 
filter, passing only those frequencies in the neighbor- 
Then, re- 


hood wo + dwo [ef. Eq. (2-13) and Fig. 5]. 


calling the phenomenon of “‘beats,’' wherein the syn 
thesis of two harmonic waves of approximately the 
same frequency gives rise to a wave oscillating with 
the mean of the two frequencies and having an ampli- 
tude envelope that fluctuates at a rate equal to the 
difference frequency, it may be inferred that the syn- 
thesis of a large number (tending to infinity) of har- 
monic waves having frequencies in the approximate 
band a + 6) and randomly distributed phases will 
give rise to a wave oscillating with frequency w) and 
having an amplitude envelope that exhibits a random 
fluctuation, the rapidity of which must be of the order 
of bw». 

The probability distribution of this envelope (as op- 
| distribution of the instantaneous re- 
sponse) may be determined by synthesizing it from the 


posed to the 


randomly oriented vectors representing the individual 
harmonic waves. This problem, known as the two- 
dimensional ‘‘random walk,” was first solved by Ray- 
leigh'® and may be regarded as the antecedent to mod- 
ern, physical applications of stochastic methods.* 
Rayleigh’s analysis gives as the probability that the 
envelope y,,(f) lies between y and y + dy 


sie y dy i 
3-4 


Ply <y,(t)<y+dy{ =P, dy =e 


rhe statistical interpretation of a sine wave having a ran 
domly varying envelope in terms of Rayleigh’s problem is due 
originally to Rice (unpublished lectures, California Institute of 


Technology, May 8, 1952); cf. also Section 3.7 of reference 1 


NDER RANDOM LOADING BY, 


This quantity P,,, known as a Rayleigh distribution 


density, is plotted in Fig. 7. 
In applying the results of Eqs. (3-1), (3-3), and (3-4 
to the essentially resonant response of a simple os- 


l 


cillator, an approximate calculation of Yo may be carried 


out along the lines of Eq. (2-13)-—viz., 
(; ) [ dw 
af a t( . 
¥ Ss 
I J 0 wW 2 fw 
[ ( | | 13° ( ) 
wo) | w 
buy 
= (})( | 3 ya) 
wy 


Yo = 


T wy Tf (a = 

= — vo [1 + 0(8*) | 3-5b 

16 F J x 
The fact that the mean square response is inversely 
proportional to 6 (rather than 6°) is a consequence of 
the fact that the fraction of the total power accepted by 
the oscillator is proportional to the product of the 
curve (cf. Fig. 


“bandwidth,” 26a), of the resonance 


»)) and its peak value (1 46°), rather than simply to the 
latter as in the case of resonance under monochromatic 
excitation. 

2-14), y 


If g(w) were given by Eq. would be given by 


V ) Ow 
vy, = fiw) dw > 3-( 
¥ “ J\e& GW, ! a-Dada 
F 70 wy 


= (F°'F,?) 4 3-6b 

yo? for Fy = VF 3-6e 

Thus y) would appear as the static displacement that 
would be produced by the root mean square loading 
Fy, = \V F*(t) in the absence of dynamic effects. Simi 
larly, s) may be defined as the static stress produced 
by Fy. These definitions will be adopted throughout 


the remainder of the analysis; however it should be re 


marked that made 


these choices of yw, Fy and sp are 
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The probability density of the maxima for essentially 
resonant loading, as defined by Eq. (3-4). 


Fic. 7 
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s=AMPLITUDE OF REPEATED 
STRESS 


log s, 
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log s-* 























log N- 
N=NUMBER OF CYCLES OF LOADING THAT 
WILL PRODUCE FATIGUE FAILURE 


Logarithmic plot of typical (s, .V) curve, as given by 
Eq. (4-2). 


Fic. 8. 


purely in the interest of convenience and in no way 
imply the approximations of Eqs. (2-14) and (3-6). 


(4) PROBABLE TIME AND EQUIVALENT STRESS FOR 
FATIGUE FAILURE 


The usual s, V curve gives the number of cycles (.V) 
of complete stress reversal (i1.e., —s to +s) of fixed 
amplitude (s) required to produce failure for a given 
structure. For this type of loading, it is sometimes 
found that there exists a minimum stress, designated 
as the endurance limit (s,), below which fatigue failure 
never occurs (i.e., VN = o if s < s,). However, the 
existence of an endurance limit depends on certain 
beneficial changes in the material produced by alter- 
nating stress at a low level (s < s,), and intermediate 
excursions to higher stress levels may counteract these 
changes; moreover, it appears that only a few materials 
(notably steel) possess well-defined endurance limits. 
Thus, it seems reasonable to ignore the possibility of an 
endurance limit when, as in the problem at hand, stress 
amplitudes are distributed over a wide range. [In 
any event, the introduction of an endurance limit has 
no important effect on the results of Section (5), pro- 
vided that the root mean square stress is large com- 
pared with this limit. | 

Experimental data usually exhibit a wide scatter on 
an s, N plot,* but within the limits of this scatter an 
adequate approximation is furnished by (cf. reference 
13) 


Nis) = (s,/s)* (4-1) 


where 5s, is the stress at which Eq. (4-1) predicts failure 
occurs in one complete cycle [although such an extrap- 
olation of the approximation, Eq. (4-1), is unwar- 
ranted, so that this interpretation of s, is essentially 
hypothetical]. On the assumption of an endurance 
limit s,, s could be replaced by (s — s,) in Eq. (4-1). 

* Indeed, it appears that N is itself a stochastic function of 
5.1% 18 
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It is customary to plot the s, .V curve on semi-log 
paper, but the approximation of Eq. (4-1) is more 
conveniently plotted as the straight line 


log s = log s; — (1/a) log N (4-2 


the slope of which is a~', as shown in Fig. 8. A typical 
value of the exponent a is 10, although it may range as 
high as 25 (reference 13). 

In order to determine the conditions for fatigue fail- 
ure under random loading, it is necessary to establish 
the ‘“‘damage’’ done by a number of stress reversals 
(7) less than the number (.V) that produces failure at 
a given level (s). The assumption 
(“‘Miner’s rule’’!*) is that the damage accumulates 


most common 
linearly, so that if a structure is subjected to n; stress 
reversals at a level s; the partial damage is ,/N(s;), and 
the cumulative damage (D) is given by 


n 
(*), x= 
i AN, 


N(s;) (4-3) 


Dy nad Z 
with fatigue failure occurring at D = 1. Assuming the 
(s, NV) law of Eq. (4-1), Eq. (4-3) becomes 


Dy = > n,(S;/5i)* (4-4 


Miner’s hypothesis, Eq. (4-3), is subject to such ob- 
jections as the fact that it ignores the temporal se- 
quence of loadings. Various alternative hypotheses 
have been offered,'® but it may be argued that presently 
available fatigue data scarcely warrant the additional 
complications insofar as the prediction of the fatigue 
life of a single sample is concerned. Moreover, as shown 
in the Appendix, the introduction of cumulative dam- 
age hypotheses more closely in accord with appro- 
priately averaged experimental results has no appre- 
ciable effect when, as in the present application, the 
stress amplitudes are continually distributed over a 
wide range. 

In formulating design criteria, it is convenient to 
introduce an equivalent stress (often called ‘‘reduced 
stress’), s,, which would produce the same fatigue 
damage as the spectrum (7, 51), (#2, S2)...after the same 
(total) number of loading. On Miner’s 
hypothesis, it follows from Eq. (4-4) that 


Sry = (> nisi" z n;|' ” 
1 


1 


cycles of 


(4-5) 


It is evident that, without postulating or otherwise 
determining an expression for cumulative damage, one 
might regard Eq. (4-5) as a working hypothesis in which 
the weighting constant could be determined either 
empirically or on the basis of some physical model 
(with or without direct reference to the slope of the 
Thus, Shanley’s theory of fatigue'* leads 
However, 


s, N curve). 
to the replacement of a by 2a in Eq. (4-5). 
presently available data are not sufficiently precise 
to sharply delimit the weighting constant, and it there- 
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STRUCT UBEL FATIGUE 


fore seems advantageous in the present analysis to write 


uf war/2a. (4-6) 
i 
and leave the choice of k open, with 2 as a (probably 
conservative upper limit. 

In the present calculation of probable fatigue life 
of a lightly damped structure under random loading,* 
the probable number of cycles of loading having an 
amplitude in the range (s, s + ds) is obtained by multi- 
plying the total number of cycles (at frequency w) 27 
by P,,(s) ds. It then follows from Eq. (4-6) that the 
probable equivalent stress is given by 


> « id l/ka 
| } gs" Pts) ds P.(s) as| (4-7) 
0 J0 


the equivalent 


Ss, = 


The (angular) frequency at which 


stress s, may be supposed to oscillate is, of course, wo. 
(5) FATIGUE FAILURE OF RESONANT STRUCTURE 


If the probability distribution P,,, for a resonant struc- 
ture is substituted in Eq. (4-7) from Eq. (3-4), the re- 
sulting equivalent stress is 


pe « . l/ka 
kat 2/2y 2y - 
5 = | sete " ds se ds (5-1 
/J 0 7 


Introducing the Gamma function (T’) 
: } gett gl ds = (5-2 
0 


Eq. (5-1) goes over to 
s, = [[(Ra/2 + 1)]! *® (QW) (0-35) 


Here and subsequently yy * represents V s°, the root 
mean square stress. 


If (since a is always large) the Gamma function is 


approximated by Stirling’s formula—viz., 
. 22+ (1/2 ‘ 
Nis+1)“e 2 V 2a, z>1 (5-4 
the equivalent stress becomes 
1/2ka) kas* : kas* ' 
Ss; = Geka)*’* > ‘ a> | 
e e 
(5-5 


It appears from the foregoing that the determination 
of the equivalent stress under random loading depends 
essentially on the determination of the root mean square 
stress produced by the same loading. If the static 
displacement yo of Sections 2 and 3 is now replaced 
by the stress (so) that would be produced by the root 
mean square force Fy under static loading conditions, 
the mean square stress produced in a resonant struc- 


ture by the random loading is given by Eq. (3-5b) as 


* The implied assumption of completely reversed loading is 
consistent with the result that the rate of fluctuation of the enve 
lope is of the order of 2 dwo, so that successive maxima and minima 


can differ only by a fraction of order 6 


U 
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7 wot (wy) 
16/ F,° 

The ratio of equivalent to static stress then is given 


by Eq. (5-5) as 


S kra wot (wo F om 
= . D-4 
So te6 F,? 


which may be regarded as the central result of the 


foregoing analysis. 


In order to form some estimate of the magnitudes 
given by Eq. (5-7), the equivalent stresses based on 
the representative power spectral density distributions 
of Eqs. (2-8) and (2-9) will be evaluated. Substituting 
these distributions, where now F? = F,2, in Eq. (5-7 


yields 


Ss, "i ka fy : ) 2 - 
— - e~ (1/2) (w/e 5-8 
So/ a (2e) 0 w) 


These results are plotted in Fig. 9. 
It is probable that Eq. (5-9) can be accepted as repre- 


sentative of random load distributions of the type 
encountered in practice. The corresponding upper 
bound for the equivalent stress is given by 

(S,/S0) max. —~ 0.51(Ra/ 5) 9-10 
For a = 10,k = 2, and 6 = 0.02 (a typical value for 


structural damping of a panel) this last result vields 
16, while for a = 0.01 (probably the 


worst possible case that could be obtained, even under 


25, k = 2, and§ = 
laboratory conditions) the ratio is 36. 
(6) JET BUFFETING OF A PANEL 


The results of the preceding sections now will be illus- 
trated by application to the specific problem of a panel 
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Fic. 9. The stress ratios given by Eqs. (5-8) and (5-9). 
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P P, (x7d)? 
LIGHTHILL'S THEORY 


w ~U /d 
eo. a 
2 2 2 i 2 2 2 
swat S?%d? x” U? exp(-Ca?*/U?) 
Fic. 10. The configuration discussed in Section (6) 


excited by aerodynamic noise (as from a jet). Prob- 
ably the most reliable procedure for determining the 
equivalent stress would be to make a direct measure- 
ment of the root mean square stress in the panel under 
the random loading of the jet (or other source) and sub- 
stitute in Eq. (5-5). If this is not possible, it might 
be possible to measure or estimate the parameters wy, 6, 
and f (wy), calculate the corresponding values of Fy and 
So, and substitute in Eq. (5-7). However, it generally 
would be difficult to estimate the power spectral density 
of the jet 7 the neighborhood of the panel, and its meas- 
urement there probably would be even more difficult 
than the direct measurement of the root mean square 
stress that it produces. Accordingly, it appears worth 
while to formulate an approximate similarity expres- 
sion for the equivalent stress on the basis of the fore- 
going results and on Lighthill's theory® of aerodynamic 
noise generation. * 

According to Lighthill's results, the density varia- 
tions arising in an aerodynamically generated noise 
field (Fig. 10) are asymptotically (¥ >> d) proportional 
to po(l'/ao)* (d/x) where l is the jet velocity (with 
respect to the surrounding medium), d a characteristic 
length (e.g., the jet diameter), . is the distance from 
the source (say the jet orifice, although directional 
effects also must be considered), py) is the mass density 
of the medium into which the sound is radiated, and 
dy is the sonic velocity in this medium. Accordingly, 
the root mean square pressure may be expected to vary 
like 
(6-1) 


V (pb — po)? ~ (U/ao)* (d/x) po 


(Lighthill's theory is not valid if 1° > ay.+) The cor- 
responding static stress in a thin plate of thickness h 
and area A then will vary like 


* It should be emphasized that, in the absence of further ex 
perimental data, the results of this section must be regarded as 
rather less reliable than the results of Sections (2)-(5 

t Preliminary results (reference 9, Part II) indicate that jet 
noise may increase much more rapidly with jet velocity when the 


latter exceeds dy. 
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A U\4 /d A 
Sy ~ \ » ie bh, ) = ™ ay { )-2 
"s 5 (.) (‘) (2) | 


The resonant frequency of a plate of thickness h, 


area A, Young’s modulus /, and specific weight w will 
vary like 


Py op " 
(gk zw 0-. 


WW ~ (h A 


while the characteristic frequency of the jet will vary 


approximately like} 
o™ U/d (6-4) 
whence 


~ (hd A) (gk wl) (6-5 


(Wy Wy] 


Substituting Eqs. (6-2) and (6-5) in Eq. (5-7) yields 


[since ayf(wo) Fy? depends only on wy w; | 


mos) GM) Le) Gd] 68 


where ® is an undetermined function. 

On the basis of both theory and experiment," it 
appears that the power spectral density of the jet has 
approximately the distribution f,(w) given by Eq. (2-9). 
Thus, if so and ww; from Eqs. (6-2) and (6-5) are posed 


in Eq. (5-9), Eq. (6-6) takes the more explicit form 


Sp C (**) ( *) ( d )() 
Po re Wy” A ‘‘h’*x Qy 
. t £e hd\? it 
exp | —Cz z (6-3 

wl A 


where C; and C, are undetermined, positive constants 
(v.i.). 

The principal shortcomings of the expressions (6-6 
and (6-7) probably are due to the omission of the 
(spacewise) phase effect in determining the loading 
of the single degree of freedom structural model and 
in the neglect of directional and Reynolds Number 
effects in the description of the jet noise. These fac- 
tors can be taken into account by introducing the 
scale factor dA in ® and C, (so that the latter is no 
longer constant) and the jet Reynolds Number and the 
angle between the jet axis and the line joining the jet 
exhaust and the panel in ®, Ci, and C.. It also should 
be noted that directional eTects become more important 
when referred to a moving reference frame, as with an 
aircraft in flight (cf. reference 9). 

Having s,, the fatigue life of the panel is given by 
N(s,) 7), where 7) is the period of the resonant mode.** 


t Experimental results vield roughly #, = al’/d (reference 9, 
Part II 

** Thus Eq. (6-7) would predict a fatigue life proportional to 
Ue’? exp (+ constant: U2 This extremely rapid decrease 
of fatigue life with jet velocity has received qualitative con 
firmation in at least one case, where the fatigue life of a panel 
decreased from approximately 100 hours to 20 min. with a nomi 


nal increase of jet power 
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STRUCTURAL FATIGUE 


In view of the great sensitivity of .V to small changes 
of s, (and, therefore, to small changes in each of the 
parameters affecting s,), it appears to be more expedient 
to base design criteria on s,. Moreover, it should be 
specifically remarked that the decrease in 7) due to 
stiffening the panel (e.g., bv increasing /;) usually would 
have a negligible effect on the fatigue life .V7\, compared 
with the large reduction of s,. Thus, while 7} is in- 
versely proportional to the panel thickness /, .V(s 


varies like ~*" exp (+ C ah*), where C is a positive 
constant. 

The damping of a vibrating panel will be due pri 
marily to radiation and structural hysteresis, but, 
insofar as 6 is expressed as a fraction of the critical 
damping and is small, the end results should not be 
too sensitive to the particular mechanism. If g, is the 
usual structural damping factor [such that the effective 
stiffness is multiplied by the factor (1 + 7g,) for har- 


monic motion |, then 


=) 25 (O-S 


6 = (1 


I 


Radiation damping is more complex, but for a panel 


it may be approximated by 


») 


6 = (2/r") (pogV A/wh) (6-9 


The effective value of 6 then should be obtained by 
adding Eqs. (6-S) and (6-9)—-viz. 


6-10 


6 => (g,/2) + (2/r° 


{ po? \ A wh 


APPENDIX 
Suppose that fatigue damage accumulates according 
to 


D = > (n,/N A-] 


where the exponent 1; depends on the amplitude of stress 
Then the probable damage per cycle re 
, the total 


reversal. 
sulting from a probability distribution P,,,(s 
integral of which is unity, is given by 


** oD 
D, = | > PGs (A-2a 
7 0 nN 
> v(s D' 1/x 
= | \y P (S§ ds (A-2b 


The damage rate, D, may be obtained by dividing Eq. 
(A-2) by the period 7). 
to be independent of s, Eq. (A-2b) vields 


In particular, if « is assumed 


OD 


Sa ite 
= ( — | : oe! alate A-3 
Ot Ty) J0 A 


Integrating this result between D = 0 and 1, the prob- 
able fatigue life is given by (using a zero superscript to 


denote that x has been assumed independent of s) 


UNDER 
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independently of the value of x. 
A more realistic assumption as to the dependence 
of x on s is furnished by recent results of Marco and 


Starkey,”’ which suggest that 
v = 2, ‘<-% (A-5a) 


r = 1, : > s, (A-5b 


where s» 1s a certain critical stress level (35,000 Ibs. per 
sq.in. for 765-T61 aluminum) at which a change in the 
Substituting Eq. (A-5) in 


type of fracture occurs. 


Eq. (A-2) and dividing by 7) then yields 


oD l (J ot] 4 ( et] 
= = = . es D (A-6) 
Or To \e \ 7, uve \ 


which may be rewritten 


OD /dot = (1/T,'”’) [A + 201 — ADD”) (A-7 


where 7,'"’ is given by Eq. (A-4), and 


[ P.,,ds [ P,,,ds 
A= : : (A-S 
‘ N / Jo N 
Integrating (A-7) then vields for the fatigue life 
7 7.0 [ dD ; 
= ( -Ya) 
* Jo A+ 2%1 — AD vate 
i 1 2-—A 
=], (1—A fi - in ( ) 
2(1 — A A 
(A-9b) 
The expression for A given by Eq. (A-S) may be 


evaluated explicitly in terms of incomplete Gamma 
functions if P,, and .V are substituted from Eqs. (3-4 
and (4-1), but, whatever the actual forms of P,,, and .V 
(except, of course, that they must be positive), A must 


lie between 0 and |. In this range Eq. (A-9b) yields 
the following values 
{ 0 0.1 02 oO 4 06 O.8 0.4 1.0 
l Ut 0 0.071 0.100 0.104 0.089 0.055 O.033 0 


Lf} 


In view of the relatively small change in fatigue life 
compared with experimental scatter) associated with 
the assumed, representative variation in exponent of 
fatigue accumulation, it is evident that Miner's hypo 
thesis (* = 1), which leads to the fatigue life 7,‘"’, may 
be regarded as adequate for the prediction of probable 
fatigue life of a structure subjected to random loading 
Further, it should be emphasized that this conclusion 
is independent of the distribution P,,, insofar as the latter 
is continuous /otherwise the differentiation (OD /On) in 
Eq. (A-2) would be invalid]. 


valid for other possible variations of the exponent x; 


This conclusion is also 


for example, if 2 is replaced by (1/2) in Eq. (A-5a), 
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the end effect is only to replace A by (1 — A) in pre- 


dicting 7p. 
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Effects of Certain Steady Motions .. . 


Continued from page 752) 
j g 


(1) The lateral oscillation becomes more lightly 
damped and the lateral and longitudinal modes of mo- 
tion are coupled at high steady turning rates. 

(2) The lateral oscillation is better damped but the 
damping of the longitudinal short-period oscillation 
becomes seriously reduced for high steady rolling rates. 
The lateral and longitudinal modes of motion are 
coupled also in steady rolls. 

(3) Damping of the lateral oscillation is reduced by 
high values of steady longitudinal deceleration. 
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STABILITY DERIVATIVES USED IN SAMPLE 
CALCULATIONS 


APPENDIX 


The stability derivatives used in the sample calcu- 
Table 2. 


values apply to a typical fighter-type airplane in hori- 


lations of this paper are listed in These 


zontal flight at three values of airspeed. All aero- 
dynamic cross-derivatives such as C,,, and C,,, are as- 
sumed equal to zero. For convenience, standard 


NACA dimensionless notation’ is used in this listing. 


TABLE 2 


Stability Derivatives Used in Sample Calculations 


Flight 
Condition V; Tt 
l 600 1.000 
2 800 0.750 
3 1,000 0.600 
Kx, = 0.1448 
Ky, = 0.4299 
€ = 1.22 Kz, = ().2686 
a = €,/2.60 Cig = 0.023 — 0.20 Cr 
CL, = 2.60 Cag = 0.060 + 0.050 C; 
Cp = 0.010 + 0.264 C,?2 Crz, = 0.40 
Cn, = —0.10 CG. = 0.035 + 0.180 C 
Cr = 1.30 Cc. = ~§000 
q 
Cn = 1.00 Cy, = 0.10 
Cmpa = Cr, = —0.20 
C.. = 0.020 0.220 Cy 
p 
Cy =0 
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An Experimental Determination of Flat Plate 


Recovery Factors 
Between lI. 


R. H. SHOULBERG,? J. A. F. 


for Mach Numbers 
90 and 3.14 


HILL,? ann M. A. RIVAS, JR.** 


Massachusetts Institute of Technology 


ABSTRACT 


An experimental investigation of temperature recovery factors 
on a flat plate was carried out in a continuous type wind tunnel 
at the Naval Supersonic Laboratory of the Massachusetts Insti- 
tute of Technology. The Mach Number range for the tests lay 
between 1.90 and 3.14, and the Reynolds Number based on the 
plate length extended to 17 X 10°. The stagnation temperature 
was nominally 110°F., while stagnation pressures were varied 
from 6 to 31 psia. The 
lows: (1) The recovery factor for the laminar boundary layer ona 
flat plate is independent of Mach Number and equal to 0.850 + 
0.007, a value about 1°; higher than the square root of the Prandtl] 
Number based on the plate surface temperature and 2°; lower 
than the square root of the Prandtl Number based on the stream 
temperature. (2) The recovery factor for the turbulent boundary 
layer on a flat plate is independent of Mach Number, but de- 
creases slightly with Reynolds Number from a value of 0.882 at 
a Reynolds Number of 5 X 10° to a value of 0.880 at a Reynolds 
Number of 10°. These values of the recovery factor are about 
1°; lower than the cube root of the Prandt] Number based on wall 


test results can be summarized as fol- 


temperature. A maximum systematic error of +0.003 is pos- 
sible, but this error would not affect the slight decrease of re- 
covery factor with Reynolds Number. (3) The beginning of 
boundary layer transition occurs on the plate at Reynolds Num- 
bers between 1.0 and 2.5 X 108. The recovery factors in the 
transition zone are higher than those of the laminar and turbulent 


regimes, but never exceed 0.90 


NOMENCLATURE 


A = area, ft.? 

( = specific heat at constant pressure, B.t.u./Ib. °F. 

h = convective heat-transfer coefficient, B.t.u./hr. ft.2 °F 
hr. ft. °F 


k = thermal conductivity of ain, B.t.u 
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l = fixed distance from leading edge measured in x direc- 
tion, ft 

\f = Mach Number 

p = static pressure, psia 

po = stagnation pressure, psia 

Pr = Prandtl Number = cpu/k 

q = rate of heat transfer, B.t.u./hr 

r = recovery factor = (Taw — 7)/(Tu 7 

Re, = Reynolds Number = pux/p 

t = excess temperature of heated nose, °F 

7 = static temperature, °F. abs 

Taw = adiabatic wall temperature, °F. abs 

7; = temperature indicated by probe, °F. abs 

7, = stagnation temperature, °F. abs 

7. = surface temperature downstream of heated length of 
plate, °F. abs 

uu = velocity parallel to plate surface, ft./hr 

x = distance from leading edge along plate surface, ft 

y = distance perpendicular to plate surface, ft 

y = ratio of specific heats 

6 = boundary-layer thickness, ft 

6* = boundary-layer displacement thickness, ft 

0 = boundary-layer momentum thickness, ft 

m = absolute viscosity, lb./hr. ft 

p = density, lb. /ft.* 

Subscripts 
( ), = free stream condition 
( )aw= an adiabatic wall condition 
Superscripts 


( )!' = boundary layer property behind normal shock 


INTRODUCTION 


Shes SUBJECT OF aerodynamic heating has received 
so much attention in the recent literature that 
further discussion of its importance seems almost un- 
necessary. Suffice it to say that wind tunnels are being 
planned with stagnation temperatures of 1,500°F. at 
high Mach Numbers, and that aircraft are being built 
whose body surfaces, unless cooled by radiation or in- 
ternal refrigeration, would reach steady-state tem- 
Missile skin tempera- 
order to 


peratures of more than 200°F. 
tures, of course, far exceed this value. In 
predict cooling requirements for the tunnel walls and 
the aircraft and missile surfaces, accurate knowledge 
of thermal effects at high speeds is necessary. Indeed, 
an acute need is currently felt for such information by 
virtually everyone interested in supersonic aerody- 
namics. 

At the present time extensive theoretical treatments 
of thermal effects in the laminar boundary layer are 
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Fic. 1. Adiabatic Flat Plate 

obtainable. The theory of turbulent boundary layers 
is much more restricted in nature, and theories of transi- 
tion are quite rare. Although many experimental pro- 
grams are now in existence, the amount of published 
experimental data relating to heat transfer at high 
speeds is greatly limited, as was pointed out by Kaye! 
in a recent survey of friction factors, recovery factors, 
and heat transfer coefficients for supersonic flow. 

The work which is described in the present report rep- 
resents the first part of a research program undertaken 
at the Naval Supersonic Laboratory of the Massachu- 
setts Institute of Technology to investigate certain 
aspects of aerodynamic heating of importance in wind 
tunnel design. The initial effort pertains to the tem- 
perature of insulated flat surfaces in high-speed flow 
or to the so-called “adiabatic wall temperature.’ This 
temperature is customarily related to the stagnation 
temperature and static temperature of the undisturbed 
(free stream) flow by means of a dimensionless recovery 
factor, r, where 

T,. — T, 
= (1) 

Tn — 1; 
Reference | showed that certain discrepancies exist in 
the recovery-factor data available in the current litera- 
ture, notably the disagreement between the values of 
the recovery factor measured for laminar boundary 
layers on cones and on flat plates, and the wide range 
of values of the recovery factor in the boundary-layer 
transition region reported on plates, cones, and cone- 
cylinders. It was with the hope of clarifying some of this 
conflicting data, as well as of extending the Mach and 
Reynolds Number ranges for which recovery-factor 
data are available, that the flat-plate recovery-factor 

tests described herein were initiated. 


EXPERIMENTAL APPARATUS 


Wind Tunnel: The tests were conducted at the 
Naval Supersonic Laboratory of the Massachusetts 
Institute of Technology, in a continuous type wind 
tunnel having an IS8- by 24-in. cross section for Mach 
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Numbers less than or equal to 2.5 and an IS- by IS-in, 
section for higher Mach Numbers. The nominal 
stagnation temperature was maintained at 110°F. 
while stagnation pressures were varied between (6 and 
31 psia. The Mach Number range covered in the tests 
lay between 1.90 and 3.14; the Reynolds numbers based 
on distance from the leading edge of the plate extended 
to 17 X 10°. Complete details of the tunnel and its 
operation are given in a descriptive manual.” 

Adiabatic Flat Plate: The design of the flat-plate 
model used in the recovery-factor tests was undertaken 
with four major objectives in view: (1) Conduction 
of heat through and along the measuring surface was 
to be minimized, the ultimate goal, of course, being the 
condition of zero heat transfer. (2) The measuring 
surface was to be a smooth, flat, continuous plane. 
(3) The model was to be mounted in such a fashion that 
its angle of attack could be rapidly altered, in order 
to set it in such a position that the strength of the bow 
wave was minimized or to change appreciably the 
Mach Number over the measuring surface for any set 
of nozzle blocks. (4) Thermocouples were to be 
mounted as close as possible to the measuring surface. 
In addition to these considerations it was necessary to 
make the plate sufficiently strong to withstand the 
maximum foreseeable air loading. All these require- 
ments were satisfied by a model with a sandwich-type 
construction, mounted in a pair of rotatable disks 
fitted into two windows on opposite sides of the tunnel 
test section. 

The model construction, shown in Fig. 1, consists of a 
'/y-in. steel plate which provides a rigid supporting sur- 
face for a '/,-in. layer of Rubatex, a multicellular plastic 
insulating material having a thermal conductivity of 
0.023 B.t.u./hr. ft. °F. and a compressive strength of 
200 psi. The Rubatex insulates a !/30-in. sheet of pol 
ished stainless steel, which provides a smooth continuous 
measuring surface, while the combination of rubatex 
and stainless steel reduces heat transfer in a direction 
parallel to the air flow to negligible proportions. The 
stainless steel sheet is fastened to the rubatex and to a 
steel frame surrounding the rubatex with a catalytically 
cured high-strength bonding material, Shell Chemical 
Corporation's Epon VI. The nose of the model consists 
of a 15° wedge of steel and Textolite. 

The plate is supported by two rotatable steel discs, 
which replace the tunnel windows ordinarily used for 
observation purposes. The discs permit a variation 
of +5° in the model's angle of attack, and three win- 
dows in the dises allow observation of the model and 
Schlieren pictures of the bow wave. A view of the 
model as seen through the open test section door ap- 
pears in Fig. 2. Further details of the constructiot 
and mounting of the plate appear in reference 3. 

Three sets of nozzle blocks were used in the tests, 
giving nominal test section Mach Numbers of 2.0, 2.5, 
and 3.0. By varying the angle of attack of the plate, 
within certain limits, surface Mach Numbers other than 
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DETERMINATION OF FLAT 


the nominal test section Mach Numbers could be ob- 

tained. 
Instrumentation: 

plate surface were made with 30 copper-constantan 


Temperature measurements on the 


thermocouples soldered on the bottom side of the stain- 
less steel sheet. Twenty-two of these couples were 
used to obtain data on the centerline of the plates; the 
other eight served to check the spanwise temperature 
distribution. Three thermocouples between the steel 
bed plate and the Rubatex were used to measure the 
temperature gradient in the Rubatex at three loca- 
tions on the plate centerline. Twelve thermocouples 
mounted on a screen in the tunnel stilling section were 
used to determine the average stagnation temperature 
of the air stream. 

The static-pressure distribution over the plate sur- 
face was measured with twenty-one 0.020-in. pressure 
taps. Thirteen taps were used to obtain data in the 
region of the plate’s centerline; eight taps, to determine 
the spanwise pressure distribution. 

Two movable boundary-layer probes were used to 
measure boundary-layer profiles. One probe, fitted 
with an external drive, was permanently mounted on 
the plate centerline with its two tips, one for impact 
pressure measurement and the other for stagnation 
temperature measurement, 2S in. from the leading edge; 
it may be seen in Fig. 2. The other probe, fitted with 
a tip for measuring impact pressure, could be inserted 
in the tunnel sidewall. The tip was located 6 in. down- 
stream of the plate’s leading edge and 3 in. from its 
centerline. 

The stagnation-temperature 


after a design of Spivack.? 


probe was fashioned 
It consists essentially of 
an iron-constantan thermocouple made of 36 gage wire, 
located in an 0.020-1n. I.D. tube which is vented to allow 
a low velocity flow past the thermocouple junction. 
After calibration in streams of known properties, the 
instrument can be used to determine local stagnation 
temperature as a function of indicated temperature, 


Mach Number, and Reynolds Number. 


DISCUSSION OF RESULTS 


A total of 3S tests, in which recovery temperatures 
were measured, were made at various Mach Numbers 
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and stagnation pressures. During 15 of these tests 
boundary-layer profiles were also measured, and during 
the final three tests a portion of the plate surface was 
artificially roughened to thicken the boundary layer 
in an effort to simulate flow at higher Reynolds Numbers 
than could be investigated with the available range of 
stagnation pressure and temperature. Since the major 
purpose of the report is to present recovery-factor data, 
the boundary-layer profiles (which were measured pri 
marily to identify laminar and turbulent boundary 
layers) are considered in Appendix A, and the roughness 
tests in Appendix B. Before discussing in detail the 
results of the recovery-factor tests, a few general com 
ments are in order concerning such matters as the tur 
bulence level of the tunnel, the dew point in the test 


section, and certain other relevant items. 


At the present time no measurement of either the 
intensity or scale of turbulence has been made anywhere 
in the tunnel. Schlieren pictures of the boundary lay 
ers on tangent ogive-cylinders with well polished sur 
faces have indicated transition at Reynolds Numbers 
of 4 & 108, 
occurred at Reynolds Numbers between 10° and 2.5 X 


The beginning of transition on the flat plate 


10°, which is in fair agreement with the summary of 


flat plate wind tunnel data given by Gazley.® 


The dew point in the tunnel was not low enough to 


eliminate the possibility of moisture condensation 
effects, and measurements in the tunnel indicated that 
condensation did occur. The principal effect of con 
densation (so far as the measurement of recovery fac 
tors is concerned) is a rise in stagnation temperature 
proportional to the latent heat of the amount of water 
stagnation tem 
The effect 


on plate surface temperature is a priort unknown since 


vapor condensing. This change in 


perature can be computed approximately. 


there exists the possibility of re-evaporation in the 
boundary layer. It was found that data obtained at 
different dew points could be correlated (with regar 1 
to recovery factor) by assuming that no re-evaporation 
occurred. 
to all the data was a change in 7 
corrections on the laminar recovery factors is shown in 
Table 1, later. The 


The only correction applied for humidity 
The effect of such 
maximum 


which is discussed 


TABLE | 


Laminar Recovery Factors 


Lowest Measured 
Value 
7, Uncorrected 


Dew Point 


Radiation and , 
Conduction 


Nose Heating Asymptotic 


Run M Pp for Dew Point Correction Correction Correction Value 
3 1.90 } 0.866 —() OOS —() 008 —() OOS 0 847 
7 1.97 6 0. 865 —() OO3 —() OO4 —() O06 0 852 

26 2 00 ( 0 872 —() O07 —0. 0038 —() O10 () 852 
13 3.43 6 0). S64 —() O02 —() OOS —() OO7 0 850 

31 2 46 6 0.874 —() O07 —() O04 —() 009 0. 854 
Q 2 56 ( 0.872 —() (9 —() OOS —) OO7 0 851 
16 2.43 6 0.876 —(). 002 0.008 —0.017 0.849 
xs 2 75 6 0.867 —() O11 -() OO] —() 012 0. 845 

35 2 9] ( 0 869 —() 008 —() O06 —(0 008 0.852 
17 3.01 Q 0.872 —() OO7 —() OOS —() O09 0. 851 
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change in the turbulent recovery factor resulting from 
such a correction was 0.008. 

Another comment concerns the general problem of 
measuring recovery temperatures on a surface which is 
not truly adiabatic. Since perfect insulators do not 
exist, one is forced to work with materials in which 
heat conduction is inevitable, especially near the nose 
of the plate where the recovery temperature differs on 
two sides of a steel wedge as a result of the different 
Mach Numbers. Radiation between the plate and the 
tunnel walls also tends to make the plate nonadiabatic. 
However, if heat transfer in directions parallel to the 
plate surface is neglected, and if the temperature gra- 
dient normal to the plate surface, the plate surface 
emissivity, the temperature and emissivity of the sur- 
face “‘seen’’ by the plate, and the convective heat- 
transfer coefficients are known, the heat transfer from 
the plate by convection and by radiation can be sepa- 


rately computed. Since 


g(convection) = hA(T, ; (2) 


7, can be determined if g, h, and 7,, are known. The 
convective heat-transfer coefficients were taken from 
Johnson and Rubesin.6 The maximum change in re- 
covery factor caused by these corrections was 0.008 
for the laminar region. In the turbulent-flow regime, 
since the maximum correction would not affect the re- 
covery factor by as much as 0.001, no corrections were 
applied. 

Some typical Mach Number distributions over the 
plate surface along the plate’s centerline are shown in 
Fig. 3. They indicate clearly that the flow is essen- 
tially one with zero-pressure gradient. The measure- 
ments of pressure and temperature in the regions out- 
side the longitudinal central zone of the plate indi- 
cated that the flow could also be considered as essen- 
tially two-dimensional. In certain of the runs tem- 
perature gradients of a few degrees per foot existed in a 
spanwise direction, but in a majority of the tests the 
spanwise variations were very small, and no consistent 
patterns of nonuniformity were noted. 


Computation of the Recovery Factor 


Assuming that the free-stream flow in the nozzle and 
test section is isentropic and that air can be treated as 
a perfect gas, it can be easily shown that 


i | 2 | (2) : I = (3) 
7 \ » foes l pi 
To y— 1 
— = | +( ) ane (4) 
T, 2 


Eq. (4) combined with the definition of the recovery 


and 


factor leads to 


: 1 te = 
— > 
Te (y — 1)M;2 


») 9 


- (¥ _— 1)M,? 
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which is the equation used for computing the recovery 


factor. 


Laminar Recovery Factors 

The recovery factor distributions along the plate 
centerline for some typical runs at a nominal Mach 
Number of 2.0 are plotted against Reynolds Number 
in Fig. 4. These values of the recovery factor have 
been computed from data corrected only for the effect 
of humidity. They also represent quite well the data 
collected at other Mach Numbers. 

As is shown in the sections of this report dealing with 
transitional recovery factors and with boundary-layer 
profiles, the initial decreasing portion of the curves 
corresponds to a laminar boundary layer. In all cases 
the recovery factor appears to decrease slowly with 
Reynolds Number until the sharp rise at a Reynolds 
Number of approximately 2 * 10°, which is associated 
with transition. It is known, however, both from the 
theoretical work of Pohlhausen’ and Crocco* and from 
test data such as those of Stine and Scherrer* and des 
Clers and Sternberg,'® that the laminar recovery factor 
is almost certainly independent of Reynolds Number. 

The nonconstant laminar recovery-factor data of the 
present tests can be explained as follows. Referring 
to Fig. 1, note that the Rubatex insulation begins 1 in. 
downstream of the leading edge of the plate. Forward 
of this point the plate construction is such as to allow a 
relatively large amount of heat transfer between the 
top and bottom of the wedge-shaped nose. Since the 
Mach Number is lower on the under side of the nose 
than on the upper side, the adiabatic wall temperature 
is some 10° to 20°F. higher on the under nose surface 
These conditions at the 


Korobkin,"! 


who shows that the entire nose wedge assumes a uni- 


than on the upper surface. 
nose are similar to those analyzed by 


form temperature appreciably higher than the adiabatic 
wall temperature corresponding to the airstream condi- 
tions on the upper surface. It is not the nose heating 
itself, however, but its influence on the boundary layer 
downstream which is of interest. The wall-tempera- 
ture excess over the adiabatic value cannot vanish 
instantly as soon as the heat transfer from the bottom 
of the plate becomes very small. The boundary-layer 
temperature profile adjusts itself over a length of plate 
which is many multiples of the heated length and only 
gradually approaches the profile which would exist 
on a plate perfectly insulated all the way from the 
leading edge. 

If it is assumed that the recovery factor does not 
depend on Reynolds Number, it follows that the ex- 
perimental data would approach a horizontal asymp- 
tote as the Reynolds Number approached infinity, if 
transition were delayed indefinitely. The asymptotic 
value would be equal to the constant value of the re- 
covery factor which would be measured on a completely 
insulated surface. If it is not assumed that the re- 
covery factor is independent of the Reynolds Number, 
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and if the variation of recovery factor with the Reyn- 
olds Number is not known, there is no way of finding 
such a variation from the present data, since it could 
not be separated from the effect of the heated nose. 
In view of the theoretical and experimental evidence 
and the difficulty of any other course, it has been assumed 
that the laminar recovery factor is independent of 
Reynolds Number. With this assumption it is pos- 
sible to calculate the asymptotic value of such data as 
appear in Fig. 5. Such a calculation may be based on 
a theoretical solution of the following idealized problem. 
A flat-plate flow is specified with a constant wall tem- 
perature (7), + ¢) for a short length downstream of the 
leading edge, 0 <x» </. Forx > / the condition of 
zero heat transfer is specified and it is required to find 
the temperature distribution in this region. The solu- 
tion was accomplished by a numerical procedure based 
on Lighthill’s'* equations for the heat transfer from a 
flat plate with nonuniform wall temperature, the details 


> 


of which are given in reference 3. A comparison with 
the experimentally observed decay of excess tempera- 
ture is shown in Fig. 5. 

The asymptotic values of the recovery factors for 
all runs with an appreciable length of laminar boundary 
layer are shown in Fig. 6, where they are plotted 
against Mach Number. These data have been cor- 
rected for conduction and radiation errors and for the 
Within the 
accuracy of the tests it is not possible to detect any 
variation of with Mach Number. 
The mean value of all experimental results is 0.850, 


effects of any water vapor condensation. 


recovery factor 
which is about one per cent higher than the square root 
of the Prandtl Number based on the adiabatic wall 
temperature and two per cent lower than the square 
root of the Prandtl Number based on the free-stream 
temperature, the Prandtl 
Number and temperature of Kaye, et al.'* The mag- 
nitude of the correction due to the heated nose is indi- 


using relation between 


Table | where the lowest measured values of 
with the 
An analysis which takes into ac- 


cated in 
recovery factor are compared calculated 
asymptotic values. 
count the inaccuracies in the measurements and the 
errors in the corrections applied to the data shows that 
the possible error in the final value of the recovery factor 


is about +0.007. 


Transition Recovery Factors 

An examination of Fig. 4+ shows that in each run, 
after the initial decreasing portion of the recovery- 
factor curve, the recovery factors rise rather rapidly to 
then 
off to an almost constant value. 
portion of the curve which centers around the peak 


some peak value, decrease, and finally level 


The rapidly changing 


value corresponds, as is shown later, to the region of 
transition from a laminar to a turbulent boundary 
layer. Every investigator who has recorded recovery 
factors measured in the transition region has noted 
similar peaks, although the reported magnitude of the 
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difference between the peak value and the values cor- 
responding to the turbulent 
0.01 to O.OS. 


highest value of the recovery factor is never more 


boundary layer have 


varied from In the present tests the 
than two per cent greater than the turbulent recovery 
factor and does not exceed 0.90. 

It is interesting to note that transition begins any- 
where from a Reynolds Number of 10° to 2.5 X= 108, 
It is possible that the turbulence level in the test sec- 
tion varies with stagnation pressure, Mach Number, 
and /or the amount of heating or cooling of the air. It 
is also possible that pressure oscillations in the stilling 
section have a role in initiating transition which be- 
comes increasingly important with decreasing stagna- 
tion pressure. 

The transition data are in good agreement with those 
of Stalder et al.,'' who record a peak transition re 
The 
only other flat-plate data extant are those of Slack,! 


covery factor of 0.89 at a Mach Number of 2.4. 


who extrapolated flat-plate heat-transfer data at a 
Mach Number of 2.4 to the case of zero heat transfer 
and found trausition recovery factors as high as 0.98. 
It might be reasoned from these results that a large 
temperature peak does exist for the flat plate, but that 
the 
This possibility has been investigated, and an 


heat conduction along the plate “‘smooths out” 
peak. 
analysis of the longitudinal heat transfer (given in de- 


tail in reference 3) indicates that the recovery factor 
at such a peak, if it existed, could have a value no 
greater than 0.906 as compared to a measured value 
of 0.900. 7 


Clers and Sternberg,'” who found no values of the re- 


The results of Stine and Scherrer’ and of des 


covery factor greater than 0.90 in the transition re 
gion, are in good agreement with those of the present 
investigation. 

In order to dispel any doubt that the velocity profile 
changed at the same location at which temperature 
measurements indicated transition, an impact pressure 
probe was placed adjacent to the plate surface and 
temperature distributions and impact pressures were 
measured simultaneously as the stagnation pressure 
was raised. The results of this test are shown in Fig. 
7 where recovery factors and the ratio of impact to free 
stream total pressures are plotted as functions of Reyn 
olds Number. The sudden change in the slope of the 
pressure-ratio curve indicates a change in the velocity 
profile which coincides closely with the recovery-factor 
rise associated with transition. 


Turbulent Recovery Factors 


Some typical recovery-factor data for the turbulent 
portion of the boundary layer appear in Figs. S and 9. 
Fig. S shows the recovery factor as a function of Reyn- 
Number at a Mach Number of 2.0; 


olds nominal 


+ This analysis clearly indicates that longitudinal heat transfer 
may safely be neglected in all but those portions of the plate very 
close to the nose or in the center of the transition region, thus 
justifying some of the assumptions made earlier 
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DETERMINATION OF FLAT 


Fig. 9, the recovery factor as a function of Mach Num- 
ber at a Reynolds Number of 5 X 10°. The latter fig- 
ure also indicates values of the cube root of the Prandtl 
Number based on the adiabatic wall temperature and 
on the local free-stream temperature. 

It is apparent that there is no recovery-factor varia- 
tion with Mach Number at a given Reynolds Number. 
At a fixed Mach Number there appears to be a slight 
decrease of recovery factor with increasing Reynolds 
Number (from 0.SS2 at Re, = 5 X 10® to 0.880 at Re, = 
107) which is also apparent in the data taken at other 
Mach Numbers, but this effect may be a consequence of 
boundary-layer changes in the transition region, and it 
high 


gradually disappear at Reynolds 


An analysis which takes into account the 


might very 
Numbers. 
inaccuracies in the measurements and the errors in the 
corrections applied to the data shows that the possible 
error in the final value of the turbulent recovery factor 
isabout +0.005. 

The data are ii good agreement with those of refer- 
ences 9, 10, 14, and 16 and in fair agreement with 
those of references 15 and 17. They also lie slightly 
below Ackermann’s'> theoretical value of (Pr 
applies strictly only to incompressible flow. 
and Maslen,'’ predict that in a compressible flow the 


, which 
Tucker 


recovery factor decreases with Mach Number at fixed 
Reynolds and Prandtl Numbers, as does Elser.*" This 
result is not confirmed by these tests. 


CONCLUSIONS 


In deriving conclusions from an examination of the 
experimental data one must remember that the varia 
tion of wall temperature with Mach Number in these 
tests was very small. This situation is common to all 
wind-tunnel test programs run at a fixed stagnation 
temperature but is very different from the situation 
occurring in flight at a given altitude. The variation 
in recovery factor is probably tied closely to the vari- 
ation of Prandtl Number at the wall. In flight the 
variation of recovery factor with Mach Number might 


be much larger than that observed here. The follow- 


ing conclusions apply directly only to flows with 
(Prj,.) essentially fixed at approximately 0.71. 
(1) For a Mach Number range of 1.90 to 3.14 the 


recovery factor for laminar boundary layers is 0.850 + 
0.007, a value about 1% higher than the square root 
of the Prandtl Number based on the adiabatic wall 
temperature and 2°; lower than the square root of the 
Prandtl Number based on the free-stream temperature. 
This value is in excellent agreement with the available 
data on cones and bodies of revolution for supersonic 
flow. The discrepancy between flat plate and body of 
revolution data noted by Kaye! and Gruenewald'® 
is now resolved. 

(2) For a Mach Number range from 1.90 to 3.14 and 
for Reynolds Numbers up to 17 X 10® the recovery 
factor for turbulent boundary layers is independent of 
Mach Number. It decreases slightly but consistently 
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with Reynolds Numbers, from a value of O.SS2 at a 
Reynolds Number of 5 X 10° to a value of O0.SSO at a 
Reynolds Number of 10°. These values of the recovery 
factor are about 1°; lower than the cube root of the 
Prandtl Number based on wall temperature. A maxi 
mum systematic error of +0.003 is possible, but this 
error would not affect the decrease of recovery factor 
with Reynolds Number. 

(3) Recovery factors in the transition zone, which 
starts at Reynolds Numbers between 1.0 X 10° and 
2.5 X 10°, are higher than those in the laminar and 


turbulent regimes but never exceed 0.90. 


Appendix A 


BOUNDARY-LAYER PROFILES 


Boundary-layer profiles were measured during sev 
eral of the tests in order to verify that laminar velocity 
profiles corresponded to (supposedly) laminar recovery 
factors and turbulent velocity profiles to turbulent 
boundary layers. Most of these tests were made in 
turbulent boundary layers, but a sufficient number of 
laminar profiles were taken to verify that a laminar 
flow did correspond to the ‘‘laminar’’ recovery factors. 

One of the few laminar velocity profiles measured 1s 
shown in Fig. 10, where it is compared with a theoreti- 
cal profile interpolated from the curves of Van Driest.*! 
Except near the plate surface the agreement is good, 
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Fic. 11. Velocity profiles in the turbulent boundary layer at J = 


3.0. 


which may be fortuitous, inasmuch as the zero height 
position of the side-wall probe was difficult to determine. 

Some turbulent velocity profiles are shown in Fig. 11. 
These profiles, which have been computed assuming a 


constant stagnation temperature across the boundary 


layer, show good agreement with a !/;-power curve. 


The measurements of stagnation temperature in the 


boundary layer yielded results of considerable interest, 
although of somewhat questionable accuracy. The 
stagnation temperature probe was calibrated in the 
tunnel free stream and the probe calibration factor, 
defined as (7; — 7)/(7, — 7), had shown small but 
consistent variations with Mach Number and stagna- 
tion pressure. In an effort to reduce these variations, 
which presumably resulted from conduction errors, a 
nylon tip was made and used in the measurements with 
the intent of calibrating it after the tests were com- 
Unfortunately it was damaged before cali- 
The stagnation temperature 


pleted. 
bration was completed. 
profiles, two of which appear in Fig. 12, were therefore 
computed assuming that the probe calibration factor 
measured in the free stream was constant. Fortu- 
nately any error introduced by the assumption dimin- 
ishes with decreasing boundary-layer Mach Number. 
The profiles, which agree quite well with those of 
Spivack,* show that the stagnation temperature in the 
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boundary layer reaches a peak value which is about 1° 
higher than the free-stream stagnation temperature. 


This peak occurs at a distance from the plate surface 


equal to about three quarters of the velocity boundary- 
layer thickness. 


Appendix B 
ROUGHNESS TESTS 


The final part of the test program consisted of an 
effort to simulate flow at high Reynolds Numbers by 
coating the upstream portion of the plate with particles 
of sand, thereby causing earlier transition and a thicker 
boundary layer than if transition were allowed to occur 
naturally. Only three tests could be made, and the 
results are inconclusive, but seem interesting enough 
to merit discussion. 

In these tests sand was cemented to part of the meas- 
uring surface with a mixture of Duco cement and ace- 
tone. Recovery factors for the unsanded length and 
velocity profiles near the trailing edge of the plate 
velocity 


measured simultaneously. The 


and 6-in. sanded lengths fit a !/; 


were then 
profiles for the 3- 
power curve quite closely, as do the turbulent velocity 
profiles measured at the same Mach Number when 
transition occurred naturally. This suggests that a 
fully developed flow of the natural transition type has 
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DETERMINATION OF FLAT 


TABLE 2 
Sand, 
Run M Po in r 6*, in. 9, in Ree 
16 3.03 30 0 0.880 0.124 0.0231 825 
57 3.01 25 6 0.879 0.224 0.0414 3100 
25 3 0.883 0.204 0.0385 2610 


58 3.01 


been reached.t The recovery factor data, shown in 


Fig. 13, suggest the same conclusion. In order to com- 
the the 


“forced’’ with those where transition occurred 


pare results of tests where transition was 


“Natu 
rally,’ the displacement and momentum thicknesses of 
the boundary layer at the aft probe position for certain 
of the runs are given in Table 2. 

Table 2 suggests that the recovery factor may reach 
a constant value (or even increase) at high Reynolds 


Number. However, the evidence is too scanty for 


more than speculation. 
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On Inelastic Thermal Stresses in Flight 
Structures 


ALFRED M. FREUDENTHAL* 


Columbia University 


SUMMAR\ 


Conventional design for thermal stresses is based on an elastic 
analysis of the stresses and their comparison with a limiting 
creep-rate or creep-strain and a time-dependent fracture-stress, 
manifestations of inelastic be 


both characteristic 


This contradiction between the premises of the stress 


which are 
havior. 
analysis and of the strength analysis necessarily leads to design 
procedures of considerable unreality. 

The effect, on the level of thermal stresses, of inelastic behavior 
of the structural material, with constant and with temperature 
dependent parameters, is investigated; the results are compared 
with those of the elastic analysis. The importance of design 
for inelastic thermal stresses is discussed 


(1) NOTATION 


a, Al = parameters 
( = constant 
c = specific heat, centigrade 
d = thickness 
D, D = plate stiffness 
r = volumetrie strain 
e = Hencky’s “strain intensity” 
E = Hooke’s modulus 
f = function 
G = shear modulus 
h = surface heat transfer coefficient 
J = tensor subscripts 
k = thermal conductivity 
K = bulk modulus 
kK, = nonlinear plate parameter 
M = bending moments 
Mri = flow moment of plate 
@) = activation energy 
R = gas constant 
s = hydrostatic stress 
s = Hencky’s “stress intensity” 
t = time (variable 
lo = time-parameter 
r = temperature difference (variable 
” Sinn = absolute temperature (variable 
Dog Was ie 
7... T,, 7 = surface, average, melting, critical, wall, 
ambient temperature 
a = variables, coordinates 
Ww = deflection 
= distance from middle plane, coordinate 
a = coefficient of linear thermal expansion 
B = surface heat transfer parameter, //k 
7 = inverse time parameter 
6; = Kronecker delta 
€, €; = strain components 
é = ratio of boundary of plastic zone to d/2 
n, 1 = coefficient of viscosity 


Presented at the Structures Session, Twenty-Second Annual 
Meeting, IAS, New York, January 25-29, 1954 
* Professor, Department of Civil Engineering 


v = time parameter, Cd?/k 
i> = coefficient of viscous traction 
= Poisson's ratio 
p = radius of curvature, density 
a, Oj = stress components 
o = yield stress 
T = relaxation time, 7/G or \/E 
¢, ¥ = functions 
w, 7 time derivative of w, 7 
mare a = first, second derivative with respect to ¥ or y 


(IT) INTRODUCTION 


y IS WELL KNOWN that aerodynamic heating of the 
boundary layer surrounding a flight-structure in 
supersonic flight The 
nonuniform temperature distribution in this structure 


produces high temperatures. 


resulting from the heating of the skin may cause 
thermal stresses of considerable severity; it has been 
shown that transient elastic thermal stresses are 


likely to attain and exceed the yield limit of the material 
after several minutes of flight. Yielding will neces- 
sarily modify the stress distribution and relieve the 
peak stresses; it will, however, also induce residual 
stresses on cooling which, under certain conditions, 
may be high enough to cause yielding of opposite sign. 
In a structural material which is perfectly elastic to 
the yield point and ideally plastic beyond it, such 
heating and cooling cycles are accompanied by rather 
severe cycles of alternating plastic straining which, 
by gradually reducing the ductility of the metal, lead 
to failure after a relatively small number of thermal 
evecles. This phenomenon of thermal has 
recently been studied under the simplest stress condi- 


cycling 


tions.” 

It is, however, generally recognized that metals at 
elevated temperatures show pronounced viscosity 
and can therefore not be considered elastic below the 
yield point; the importance of this deviation from 
elasticity increases rapidly with increasing temperature. 
Such recognition is expressed in the specification of a 
limiting creep rate or a limiting amount of creep 
deformation as the principal design criterion for high 
At 
the same time, however, the stresses which presumably 
cause such the 
are computed on the basis of elastic analysis. 


temperature equipment such as heat exchangers. 


stresses, 
While 


such contradiction between the basis of stress analysis 


creep, including thermal 


and of design necessarily leads to discrepancies in 


high temperature design of predominantly —load- 
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INELASTIC 


arrying elements, it results in a completely unrealistic 
subject to temperature 
sradients or variations. This is due to the fact that 
the inelastic response of the structural material only 


lesign of elements severe 


fects the distribution of the load stresses, the resultant 
‘ which is independent of such response, while it is 
the intensity of the self-equilibrating thermal stresses 
which is itself primarily determined by this response. 
In the analysis of the stability of columns at elevated 
temperatures, as well as in the analysis of load stresses 
in beams and rotating discs the effect of creep has been 
investigated for various forms of the creep law; in 
thermal stress design, elastic analysis is still the rule. 
Rational design for thermal stresses, however, should 
be based on a stress analysis consistent with the 
behavior of real materials and with the design as- 
sumptions concerning their deformation and strength 
it elevated temperatures. * 


III) Types oF INELASTIC BEHAVIOR 


The deformational response of structural metals 
it elevated temperatures is a combination of elastic, 
The three 


mechanical constants associated with these components 


juasiviscous and plastic components.’ 


{ deformation—the elastic modulus, the coefficient 


i} viscosity, and the yield stress--are temperature 
lependent, and so are the thermal constants which 
letermine the level of the 
oeflicients of linear expansion and of thermal con- 


thermal stresses: the 
luctivity. 
Of the 
variation of the 
important, particularly in relation to the very rapid 


mechanical constants the temperature 


elastic modulus is relatively un- 
change with temperature of the coefficient of viscosity. 
It is also quite likely that the apparent, rather slow 
reduction of the modulus of elasticity with increasing 
temperature is, at least partly, due to the difficulty to 
bserve the true elastic response of a material in the 
presence of that, therefore, the 
modulus is, in fact, even less temperature-dependent 


creep and elastic 
than suggested by test results obtained in conventional 
experiments. 

The coefficient of viscosity is the most temperature- 


sensitive of the mechanical constants. The relation 


1 = mexp (Q/RT*)/exp (Q/RT,) = 


wit 
nm exp | O/ RT, mo Ge | 


J 


In (n No) = 


is generally used® to describe the change of this coeffi- 


jcient from its value m at the absolute (Kelvin) tempera- 


ture 7) to its value 7 at the temperature 7” = 7) + T. 
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For temperature changes that are small with respect 
to 7, Eq. 
In (n/m) = O RT, [#7 7, 4 


(la) can be approximated by 
T/T,)?] = 
#QO/RT,(T/T, Ib 


For metals the ratio (QR) (Q in cal. per mole and 
R = 2 cal. per mole,°C. 
tude 10* to 10°. 


The change of yield stress with temperature varies 


is of the order of magni 


with the material. It is very slow for high-temperature 
alloys such as Inconel, moderate for titanium alloys, 
but rather rapid for the more conventional structural 
materials such as 75S-T6, SAE-4340 steel or stainless 
steel W. 


the observed change® can be roughly approximated is 


A relatively simple expression by which 


of the form 
Co = Boll 2T/T»)?] 2 


where & is the yield stress at the reference temperature 
7 = 0, and 7°, the melting temperature. 
The thermal 


increases 


coefficient of expression generally 


with temperature; the relation can be 


roughly represented by the function 
a= a+ AaVT/T, 3 
where a is the coefficient at the reference temperature 


T = 


0 and Aa the total range of variation between 


7 = 0 and the melting point. 

The thermal conductivity decreases for certain 
metals, while it increases for others; thus low-carbon 
steel at 1450°F. shows a decrease to one-half of its 
room temperature conductivity, while Inconel X 
doubles its conductivity over the same range; the 


variations may be assumed as roughly linear 

The use, in thermal stress analysis, of the initial, 
constant coeflicient of expansion results in an under 
and in an 


estimate of the stresses due to heating, 


overestimate of the stresses due to cooling In the 


case of the thermal conductivity decreasing with 
temperature, the use of the initial conductivity coeffi 
cient overestimates the rate of heat flow and might 
therefore underestimate the early transient stresses; 
the opposite trend in the variation of the conductivity 
will produce transient thermal stresses smaller than 
those computed on the basis of the initial coefficient of 
conductivity. While the effect, on the 


stresses, of the temperature-sensitivity of both thermal 


thermal 


coefficients may thus be significant, the consideration 


of this variation in the heat transfer and thermal 
stress equations makes their solution rather difficult. 
The use of intermediate constant values will improve 
the accuracy of the approximate solutions. 

The simplest idealization of the mechanical behavior 
of structural materials at elevated temperatures is the 
with 


assumption of viscoelastic 


to the yield limit and ideal plastic or visco-plastic 


response creep up 


response above this limit. The viscous component of 


the deformation is usually nonlinear; this complicates 








JOURNAL OF THE 
the solution of even simple visco-elastic stress problems 
to such an extent that it appears expedient, particularly 
in a general discussion of the effect of creep in design, 
This 


response can be assumed to roughly represent the 


to consider first a linear visco-elastic response. 


behavior of metals either at relatively low stresses or 
at very high temperatures; at stresses in the vicinity 
of the yield limit and within the range of moderately 
elevated temperatures the nonlinearity of the viscous 
response must be introduced. 

The linear visco-elastic stress-strain relations for a 
material with creep (Maxwell body) have the general 


(5 l ) h 
i (o — Fe 
of T 


where the relaxation time r = 7 G. 
Nonlinear viscous steady-state creep equations of 


form 


2G (e€; — ¢6; = 


different forms have been proposed; the simplest form 


is that of a power function 


where s, and e, denote Hencky’s stress and_ strain 
intensities which are proportional to the square root 
of the second invariant of the deviators of stress /» and 
strain Jo, respectively’ 


and (6) 


A, and » > | are constants because of the desired 
symmetry of the stress-strai:: rate relation » should be 
an odd integer. 

For stresses exceeding the yield point the yield 
condition 

Sy = 0% (7) 

determines the state of stress in the ideal plastic region 
of the material. 

In the linear the 
relation between the components of the stress deviator 


visco-plastic material general 


and the strain-rate deviator 


atv d= (1S Norn — 8 
2n| (€;; — e6;;) ={1—- (oi; — S0:;) (Ss, > oo) (8) 
Ol 
For s, < o the visco-plastic strain rate vanishes; the 
behavior is then governed by the visco-elastic relation 
(4), the simplifying assumption for establishing Eq. 
(S) being that the visco-elastic strain rates are neg- 
ligibly small in relation to the visco-plastic strain 
rates. If this assumption cannot be made, Eqs. (4) 
and (8) must be combined for s, > op». 

For all materials the relation between e and s has 
the form e = s/3K, 

For the discussion of thermal stresses in plates the 
above equations of inelastic response can be simplified 
Ts, = @; 0 and 


by the assumptions o,, = i = 
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therefore ¢€,, = € y = € 


/ 


Under this condition the 


following relations are obtained 
Linear visco-elastic material 


O (5 ) a 2 Oo O 
e=> > om + Z 
Or Oo T/ 6G QA Of Oo 
vo | 


1 — a O a 
‘3 — (al fa 
1+ vol 37T/2G Or 


Nonlinear viscous material 
R) | J ii ao ; 
e = -_(a/X,) + > t (al a 
Ol 2 QA Of Of 
Nonlinear visco-elastic material 
—py 1 Oa | O 


o | 
+ v2G or 2 or 


¢ = (al ob 
Or I 
Ideal plastic material beyond the yield point 
G¢ = Goa (7a 
Linear visco-plastic material beyond the yield point 
O | 2 0c O 


€ = (C0 — Go ia ~ i 
9K of Or 


) al (Sa) 
ol On, 


Linear visco-elastic visco-plastic material beyond the 


vield point 

re) (| 
€ = 

Ol 


The above equations are 


— vo l\o l Oo aa 
ia — (of — Go) +f (aT) 
+ voOt 37/2G_ 6m Of 
(Nb) 


somewhat simplified if 


incompressibility of the material can be assumed 
[A — o and vy > (1/2)]. 
(IV) SURFACE STRESSES IN PLATES SUBJECT TO 


RAPIDLY CHANGING AMBIENT TEMPERATURES 


(a) Instantaneous Change of Temperature(Thermal Shock) 


Considering an infinite slab of thickness d or a curved 
slab the radius of curvature of which is large in relation 
to its thickness, with one face exposed to a rapidly 
changing ambient temperature, the other face insulated 
or (O07 /0z) = 
is one-dimensional and equivalent to that of a plate 


0, the temperature-distribution problem 


of thickness 2d exposed at both faces to symmetrically 
If the plate is not 
restrained along the edges, the thermal stresses at 


changing ambient temperature. 


any distance z from the center plane will be the differ- 
ence between the stresses suppressing the local thermal 
expansion depending on 7(z), and the uniformly dis- 
tributed stresses due to the forces at infinity re-estab- 
lishing the force-free edges. For the linear elastic 
or visco-elastic plate with temperature-independent 
parameters, in which the thermal stress intensity is 
directly proportional to the temperature, this stress 
difference is proportional to the difference (7, — T) 





betw: 
temp 
must 

tract 

stress 
force 
maxi 
incre: 
will 1 
varia 
shift 

plate 

In 

surfa 
distri 
know 
probl 


series 


is we 
nenti 
paral 
equal 
heat 
To 
inelas 
the 1 
must 
howe 
nonli 
portic 
in th 
funct 
8, obi 
the s 
stress 
dimet 
from 
relate 
to th 
spect: 
funct 
easily 
of the 
analy 
On 
Eq. | 
funct 


where 
pre rd 
well; 

relati 


on. the 


(Sb 


\ia 


point 


(Sa) 


d the 


/*) 
(Sb) 


od if 
ined 


10Ck) 


rved 
ition 
idly 
ated 
lem 
slate 
ally 
not 
s at 
ffer- 
‘mal 
dis- 
tab- 
istic 
lent 





INELASTIC 


between the local temperature 7(z) and the average 
temperature 7,. For nonlinear behavior the stresses 
must be computed by setting (de/d/) = O, and sub- 


tracting, from the resulting stress o7(z) the uniform 


stress om = (1/2d) So aordz, produced by the 
jorce which re-establishes the free boundary. The 
maximum stress ¢ = om — or (compression for 


increasing, tension for decreasing ambient temperature) 
will usually be at the surface. Consideration of the 
variation with temperature of a@ and k may, however, 
shift the maximum stress toward the interior of the 
plate. 

In order to compute the surface stress not only the 
surface temperature, but the full temperature 
distribution over the thickness of the plate must be 
transient heat-transfer 


also 


known. The solutions of the 
problem giving the distribution in terms of an infinite 


series 
(T/T) = 1 — Yogxlt) -¥il2 (9) 


is well known; the time functions are negative expo- 
nentials ¢,(¢) = *t/3), 3) = Cd*/k is a time 
parameter and a, are the roots of the transcendental 
dh k is the surface 


exp | —< ls 


equation a,-tan a, = 8, where 8 


heat transfer parameter. 
To compute the stresses or and oy) in nonlinear 
inelastic materials by combining Eq. (9) with any of 


Eqs. (4) to (S), numerical methods 


the nonlinear 


must be applied. In a first approximation it may, 
however, be assumed that in both linear and slightly 


nonlinear materials the stress ¢ = aw 


or iS pro- 


(T, — T 


difference, as a 


portional to the temperature difference 
in this therefore, the 
function of (¢/3) for different values of the parameter 


case, critical 


8, obtained from a numerical evaluation of Eq. (9) for 


the surface at which 7 = 7,, determines the surface 
stress. If the peak values of (7), — 7) and the 
dimensionless time of their occurrence (f)/ 3), taken 


from plotted (7, — 7.) versus (¢/d3) functions, can be 
related to the parameter by fitting simple functions 
to the points [((7, — 7\)mar, 8] and [(to'd), 8] re- 
spectively, and if the general shape of the plotted 
functions can be reproduced by a relatively simple, 
easily integrable function, the effect of the inelasticity 
of the material on the thermal stresses can be roughly 
analyzed. 

On the basis of published numerical evaluations of 
Eq. (9) for various parameters 6° it appears that a 
function of the type 
(10 


(y/y)" = x-exp (1 — x) 


where y = (JT, — T7;)/To, Yo = Vmars X = t/t re- 
produces the shape of the theoretical functions fairly 
well; by curve fitting the following approximate 


relations, valid for 5 < 8 < 25 have been obtained: 


THERMAL STRESSES IN 


FLIGHT STRUCTURES 77d 
0.1S/V B; 


n= 4 1] 


vw = 1/11 + V59/B)s ty /d = 


Introducing (11) into (10) the equation for the transient 


change y(t/d 


— A(p)(t d “exp |—a(8)t/d] 12) 


where 


A = 1.96 YB (1+ V5 8) and a 


1.39 vWf B (5 < B < 25 13 
The effect of the inelasticity of the material on the 


shock 


introducing the difference 


thermal stresses under conditions of thermal 


can now be discussed by 
between the rate of deformation due to 7 and the rate 


of deformation due to the negative force restoring the 


: (a) 
free edge and roughly proportional to 7 - 
Of 
. I) < ; T, 7.) 4 he inel 
(a = (J, — 7.) into the inelasticity equations. 
oy or | | e) Ss 


Thus, for instance, for the simplest case of linear visco 
elastic response the differential equations for the stress 


at the surface os 


o i1+71  t+v0 
_ f Cs = PGa 7 Ts 14 
Ol l— vor | vol 


Assuming a temperature-independent relaxation time 
37r(1 — v)/(1 + v) and introducing (7 7's 


according to Eq. (12), the solution of (14) can be 


written in the form of a series 


l+p ‘te Fe 1/43 
Ae~*? (t/9)'41 1 Y 
l_—~y» I\T 


os = —2GaT, 


where (5) elastic 18 the expression in front of the bracket 
in Eq. (15). For the maximum temperature difference 


at ¢ = ¢, the surface stress 


j if 1 / 4b, 
r.) elastic a l 7 =e l T 
{ OT} Q\ 7, 
| (= \) l to |! 
; — _ ra 2 = ] ain 
9.13 T} 9.16.04 is |S 


The relative surface stress intensity o,(1 — » 
2GTya(1 + v) according to Eq. (15) has been evaluated 
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as a function of (t/3) for the ratios (3/7,;) = 25, 10 
and 0 (elastic case) and the surface heat transfer 
parameter 8 = 5, and presented in Fig. 1. The dia- 


grams show that for a sudden increase 7) of the ambient 
temperature the visco-elastic response of the plate 
causes a reduction of the maximum elastic stress 
intensity followed by a stress reversal prior to the 
complete stress relaxation. The sharp drop from the 
maximum stress intensity of the visco-elastic curves 
also suggests that even for delay times f) that are short 
in relation to the relaxation time of the material, so 
that almost the full elastic stresses can be built up, 
these stresses are only of short duration, decrease at 
a much faster rate than the temperature difference 
(T, — Ts) itself, and therefore disappear before the 
temperature gradients disappear. Since the maximum 
stress is of short duration, the critical design strength 
values of the material are the very short-time strength 
values, representing the stresses that can be safely 
sustained during a few minutes or even fractions of 
minutes; these stresses are significantly higher than 
the values which are usually considered to represent 
short-time strength. 

When the response of the material is pronouncedly 
nonlinear as, for instance, in the case of an elastic 
plastic plate, a plastic region will spread from the 
surface as soon as a, => o» and extend as far inward as 
the elastic stress o > oy. Subsequent reduction of the 
difference (T, — T7,) will produce elastic stresses of 
opposite sign proportional to this reduction, and thus 
residual stresses orgs = 9% Selastic: In the course 
of a rapid change of ambient temperature of such 
magnitude that the elastic stresses at the peak difference 
OF T,) exceed the yield point over a finite depth 
from the surface, the subsequent reduction of (7, — 7,) 
due to heat flow [Eq. (12)] may produce severe elastic 
stresses of opposite sign and therefore substantial 
residual stresses, even without reversal of the ambient 
temperature. If the residual stresses are not relaxed, 
increase the elastic stresses resulting from 
and at the difference 


they wi! 
such reversal before 
[—(7, ts) 
Whenever the yield stress is attained in a visco- 
elastic-plastic plate the stress will remain there only 
During this time the visco-elastic 


peak 


for a short time. 
stress distribution 
are built up as soon as the stresses decrease below the 
however, 


is changed and residual stresses 


vield point; these residual stresses are, 


eliminated by relaxation. 
(b) Aerodynamic Heating 


Conditions of aerodynamic heating are much less 
severe than those of thermal shock. It has been shown’ 
that under certain simplifying assumptions the mean 
temperature increase of the skin is governed by the 


relation 


Tw(t) = Twoll — e~*) (16) 
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if the initial temperature of the structure is taken as 
zero, Where Two is the adiabatic wall temperature rise 
h/Cpd is an inverse time parameter; in 
effect of heat 


and y = 


this relation the radiation is neg- 
lected. 
that the 


transient temperature difference between surface and 


Assuming, as a rough approximation, 
center of an idealized element of a flight structure 
produced by a temperature increase d7w of the skin 
can be described by Eq. (12), the total temperature 


difference is obtained by integration 
(T, — Ts)/Two = —Aryfi(t Bo) fe STOIC) ae (97 


Introducing this expression into Eq. (14), the differen- 
tial equation for the visco elastic stress at the surface is 


obtained 


(; | Gal we 
oo os = —2G 
Of T 1 : wes —— 


( IS 
the solution of which is 
rT l+» P 1 } ) t 
os = —Gal we 1yt(t/d)'/e |! 1 
ans i 
| ££ | 
c(t/d) + c- (td — — o°(t/3)*® + 
0) 9.13 9.13.14 
(19 
where c = (a + yd — B/7; The maximum surface 
stress is obtained by solving the equation (dos/dt) = 0; 


neglecting terms in (f 3) of higher than first order, 
the time ¢) at which this stress is reached is, in first 
approximation, (f/3) = 9/4(a + yi + 8/7). Hence 
chy /F = 91 — «)/4Q1 + k 


T} —> co 


where «x = B/r7,(a + yd). 


For elastic response and therefore x > 0; 
the visco-elastic stress reduction depends on the value 
of 0 < x < 1, in which the effects of thermal coefficients 
of the material, surface heat transfer resistance, rate 
of heating and relaxation time are combined. For 
x = 0.5 the maximum visco-elastic surface stress is 
about 0.42 of the elastic stress. 
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\V) SURFACE STRESSES IN PLATES SUBJECT TO LINEAR 
TEMPERATURE GRADIENT 


If the two surfaces of a flat plate or of a plate the 
curvature of which is small in relation to its thickness 
d are kept at temperatures 7, and 7; a steady state 
of heat flow is established after a certain time, which 
is associated with a linear temperature gradient. 
The difference between the change of length of the 
external fibers and the middle surface produces a 


spherical curvature of the unrestrained plate 
lja = ei, — T,)/d = el /d (20 
If the edges of the flat plate are perfectly restrained 
against rotation, or if the plate is part of a circular 
cylinder or of a sphere of large radius, a constant 
bending moment 1J/, = J, = J is induced which 
restores the initial curvature. 
Since for the elastic plate the relation between 
deflection w and bending moment 


Wer = Wy, = —M/D(1 + v) = —l/p 


where D = (21 


the temperature moment 

M = alTD(1 + v)/d (22 
and the maximum elastic fiber stresses 
takT/2(1 — pv (23 


c= 


For the linear visco-elastic plate according to Eq. (4a 


l1— vo i: o . 
Pr =v, = —- ~ 6\1/Gd* = 
1+ vO OT 
d , 4 
om (| p (24) 
dt 
since at the surface €,, = €,, = *dw,,/2 = ¥dw,,/2 


and o = 6M/d’. 
obtained from the differential equation 


|— vo l = , 
( + Jar = —aTlD,/d 
1+ vot 3r 


where D 


Hence the temperature moment is 


= Gd*/6 (25 


For constant 7 and 7 = 0 the solution of this equa- 
tion describes the relaxation of an initially induced 
temperature-moment .l/), which proceeds exponentially 
with the relaxation time 7; = 37(1 — v)/(1 + pv). If 
the temperature difference 7 increases slowly so that 
the temperature gradient can be considered to remain 
invariably linear, the temperature moment is obtained 
as a solution of Eq. (25) for a given time function 7). 
For 7(t) according to Eq. (16 

M(t) = aTD,ny(e ¥ — @ d (1 — yr) (26 


M=0O0ati = 
a maximum, which is attained at f = 7, In 7:7 


The moment thus increases from 


(ny-—1 


THERMAL STRESSES IN 


FLIGHT STRUCTURES Gia 


This maximum is considerably smaller than the elastic 
moment if 7; and 1 y are of the same order of magni 
tude; for 1 y = 0.57, the visco-elastic temperature 
moment does not attain one-half of the elastic moment 
For a linear increase of the temperature difference 

T = Tot /t 
Mit = al,D,7;(1 é f 27 


If, from a time ¢ = ¢, this difference remains constant, 
the moment J/(f,) relaxes exponentially with relaxation 
decrease of the 


time 7,. Therefore, a 


difference, after a time which is sufficiently long to 


temperature 


permit a substantial relaxation of the temperature 


moment caused by the preceding increase of the 
temperature difference, will produce thermal moments 
of opposite direction, as this increase is equivalent 
to an increase of negative sign with respect to the new 
stress-relaxed condition. Hence, in visco-elastic plates 


sufficiently long cycles of temperature difference 
between zero and a maximum will be accompanied 
by almost complete moment- and _ stress-reversals, 
leading to fatigue failures. 
Consideration of the temperature dependence of 
7, would require the introduction of a variable viscosity 
over the thickness of the plate according to Eq. (la 
The resulting stress distribution is nonlinear, since the 
a/MT) 


is obtained 


viscous strain rate € = d/dt. (z/p) — Sw, 
where f(z) 


and therefore ¢ = —w,,2-f(z 
by introducing a linear temperature distribution into 
Eq. (1b). With 


the nonlinearity becomes more pronounced; the neutral 


increasing temperature differences 


axis moves toward the colder surface in which the 
fiber stress tends to increase faster than in the hotter 
surface. 

For an incompressible nonlinear viscous material 
according to Eq. (5a) the differential equation for the 
thermal moment is obtained by expressing « and o in 


Eq. (5a) in terms of p and / and introducing p ac- 


cording to Eq. (20). Hence 
_,dfi - . 
Ww" = -K = —alk, /d (2S) 
dt\p 
where 
K, = nd-"*"-d,/2(2n + 1 (29) 
Therefore for the nonlinear visco-elastic plate ac- 
cording to Eq. (5b) the differential equation of the 
temperature moment .\/ 
1 — vO ae = 
(M/D,) + (M/K,)" = —aT/d (380) 
1+ vO 
or, with c, = ([(1 + »)/(1 — v)|D,/A 
j l+, - . 
M+eM" = —- al D (50a) 
l—~yp 
For a stationary temperature difference 7 const. 


and 7) = O the relaxation of an initially induced 


temperature moment ./, thus proceeds more slowly 
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than in the linear material, according to the equation 


l 


(M/M,)"~* = ee 
1+ (nm — 1)My" 


(ie > 1) (I) 


Cyl 
For a linear increase of the temperature difference 
T = Tot/to the equation 


| 


a es — 
M+e,M" = - aly)D,/t) = const. (32) 

‘ne 
is directly integrable; for less simple temperature 


histories numerical integration using finite differences 
is required. 

In an elastic-plastic plate the yield limit is reached 
in the surface fibers when oa, according to Eq. (23), 
the critical temperature difference 7), 


reaches ao; 


therefore 


Toe = 200(1 — v)/ak (3: 

If the temperature dependence of the yield limit 
according to Eq. (2) is considered, the critical tempera- 
ture 7, the temperature 7), is 


obtained from the equation 


rr 


in terms of critical 


Tcl l _— t( i . )>] < 


For aluminum with o/E = 2 X 10%, a = 10-5 
per °C. and v = 0.35, the critical temperature difference 
Toe = 130°C.; with 7,,, = 660°C. the critical tempera- 
ture difference corrected for temperature-dependent 
vield limit 7, = 113°C. 

For 7 > 7. the plate under a constant moment is 
t£(d/2) from 


in the elastic state over a depth s = 
the middle plane, and in the plastic state between 
gs = +£(d/2) and zs = +d/2. Since the curvature of 
the plate is constant and determined by the maximum 
strain € in the elastic core 

2o0(1 — v)/Héd = 


aT./td (34) 


l ; = es yy = 2€ td ani 


and since the relation between the bending moment 
Mand the depth of the elastic core in an elastic- 
plastic plate" 


a (335) 


M = M rz [1 — (1/3) 


Str 


where the ‘‘flow-moment” 1/7, = ovd*/4 denotes the 
theoretical ultimate plastic carrying capacity of the 


plate, the relation between & and 7° obtained from 
Eqs. (20) and (34) 
Ee 7/1 (36) 


determines the relation between temperature difference 
and temperature moment in the elastic-plastic plate: 


l 2 
(7 ,/")* (5 Mv ~<Méz Mr) 


3 


M Mrz =]1-—- 
(37) 


A reduction of the temperature difference is accom- 
panied by a reversed elastic moment according to 
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Eq. (22), with a temperature difference (—7 Hence 
for 7 > 7), a residual bending moment JV/pps will 
remain in the plate after complete removal of the 


temperature difference 


i 2 _ 
Mrs = M rz l = a (zy. 1 y= — ms 1 1 as 
») , 
The residual fiber stress 
SrEs = oo(l — T/T, 39) 
For T > 2Ty,, Mears > (5/12)Mrz, and ores > 
hence the removal of the temperature difference is 
accompanied by plastic deformation in the opposite 
direction. Repetition of temperature cycles in excess Fc 
of 27, will therefore gradually lead to failure by 
“alternating plasticity” or ‘‘thermal fatigue. N 
Ye 


When loads and temperature differences act jointly 
on an elastic-plastic plate the moments due to loads 
and to temperature-differences are not independent. 
Even in the simplest case of a constant load moment 
Jp and a constant temperature difference 7 producing 
M, in the 


direction of the load-moment, the dependence of 1/ 


a constant temperature moment acting 


on Jp cannot be easily evaluated. It is obvious that r 


if, in this simple case, 1/p > °/;.\/¢,, the temperature " 


moment for the same difference 7 will decrease for 


increasing load moment if the load moment is applied 


For less simple cases, particularly those of 


variable load 
different portions of the plate a highly complex inter 


first. 


moments with diferent signs over ct 


action must be expected in the elastic-plastic as well 
as in the nonlinear visco-elastic plate. Only in the 
elastic and linear visco-elastic plate with temperature 
superposition of 


independent is simple 


load and temperature effects permissible. 
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of Technology, Cambridge, Mass 
ine 4, 1954 
HE TITLE PROBLEM has received considerable attention re 
cently ind is of considerable practical and theoretical in 
terest Principal available results are well summarized by Tri 


bus and Klein 
forth by 
Klein 


which is restricted to the case of almost constant 7, the shear dis 


ind the significance of variable surface tempera 


tures is set Schlichting.2 Among their new results 


rribus and have generalized earlier work by Leveque 


tribution It is shown herein that this restriction is not essential 


since a simple coordinate transformation reduces the problem 


for arbitrary 7 to that treated by Leveque 


As noted by Tribus and Klein the energy equation (we assume 
constant fluid properties, incompressibie flow throughout 
p ¢ u 7 t+ov7,) =k Ty, + wl . l 


is linear so that, given any particular solution, the problem 1s to 
solve 
pC,(uT, +vT,) =k 2 


subject to appropriate boundary conditions. A typical set ts 


T(x, 0) = G(x given) | 
T(x = () } 
7(0, y) = 0 \ 


When it is possible to approximate u by 


T( AX 
“= bd } 
v1 
the continuity equation vields 
: re : 
v,=- y? ) 
9 : 
- 4b 


) 


If + can be taken as constant, i.e., 7’ is negligible, Eq. (2) becomes 


This problem has been solved by Leveque for a step function in 
temperature and Tribus and Klein have considered arbitrary sur 
face distributions of temperature and heat flux. The restriction 


ou 7 is often too severe and Eq (5) must be considered 


Retaining 7, Eq. (2) becomes 
pf , = 
ryl, — r y27 = k 7 (7 


Equations of this type occur in other heat transfer problems 
where a reduction to the form of Eq. (6) is made 

* This work was done as part of a program on heat transfer in the Naval 
A F33(616)-2042 with WADC 


Supersonic Laboratory under Contract 


ronaulica Clence 


. are presented in th 


and discu 
/ j >} j / 
special department Publication com pleted 


The Editorial Committee does not hold 


y lhytheco pondent 
Let 
° 
yV or wf Vado 
J0 
Phen Eq. (7) becomes 
p Cpnl k 1 8) 


Since in Eq. (6) 7 1s considered constant, it is equivalent to Eq 


S The form of the boundary conditions is not affected As 


i result the solution obtained in reference 1 can be used with arbi 


trary 7 (hence, pressure gradient) as well as temperature or heat 


flux at the wall 
i different way by Light 


Similar results have been derived in 


hil who notes that the approximation implied by Eq. (4) cor 


responds to flow of high Prandtl Number fluids 
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Hea I 


a body which suddenly starts moving in 


aND LoMAX! have pointed out that the pressures on 
a compressible 
fluid depend, at the first instant, only upon the local normal v« 
locities of the surface 

Other situations in which the gas suddenly is made aware of the 
existence of an obstruction, as for certain problems in hypersonic 
flow, have also been solved by considering locally planar waves 
emanating from the body.’ The observation made in refer 
ence 1 can also be applied to rapidly oscillating bodies, consider 
ing the frequency of oscillation as the primary variable rather 
than the time-constant of the motion 

In this note, it is shown that such a simplification of the prob 
lem of finding the pressures on oscillating wings can be deduced 
small disturbance considerations, and criteria 


directly from 


for the applicability of the simplified theory can be obtaine | 


in 


the process. The general expression for the pressures on an oscil- 
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TABLE | 
Exact and Approximate Force and Moment 
Functions for the Oscillating Supersonic Airfoil (M = 2.5) 


le L L L L 

/ Ly 2 2 3 3, ly mah 
0 0 0 0 0 0 0 0 

0.119 0.0003 0 0.0473 0.0476 0.0057 0.0057 0.0476 0.0476 

0.238 0.0008 O 0.0950 0.0952 0.0227 0.0226 0.0952 0.0952 


0908 0.0906 0.1886 0.1904 


oO 


0.476 -0.0100 0 0.1903 0.1904 
0985 0.0984 0.1963 0.1984 


oO 


0.496 -0.0108 0 0.1967 0.1984 
0.992 0.0241 0 0.3755 0.3968 0.4019 0.3936 0.3742 0.3968 
1.253 0.0411 0 0.4908 0.5012 0.6526 0.6280 0.4638 0.5012 
1.701 0.0580 O C.6938 0.6804 1.2250 1.1574 0.6178 0.6804 
1.984 -0.0641 0 0.8224 0.7936 1.6793 1.5747, 0.7160 0.7936 


2.480 0.0706 O 1.0458 0.9920 2.6441 2.4601 0.8886 0.9920 





ik : 
M) M M M2 M3 M3 M, M, 





0.119 0.0006 O 0.0478 0.0476 0.0057 0.0057 0.0635 0.0635 
0.238 0.0016 O 0.0953 0.0952 0.0227 0.0226 0.1273 0.1269 
0.476 -0.0200 0 0.1978 0.1904 0.0898 0.0906 0.2537 0.2539 
0.496 -0.0225 0 0.2025 0.1984 0.0972 0.0984 0.2643 0.2645 
0.992 9.0188 O 0.3530 0.3968 0.3904 0.3936 0.5043 0.52 
1.253 0.0447 0 0.4663 0.5012 0.6379 0.6280 0.6235 0.6683 


1.701 0.710 0 0.6712 0.6804 1.2074 1.1574 0.8281 0.9072 


1.984 0.0806 0 0.8016 0.7936 1.6607 1.5745 0.9585 1.058 
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2.480 0.0908 0 _ 1.0281 0.9920 





lating wing is developed, and the results of sample calculations 
are compared to results of more exact theory. Lin, Reissner, 
Tsien’? have laid the systematic foundations for treatment 


Extending their consider- 


and 
of unsteady flow over slender bodies 
ations to flow over wings of finite aspect ratio, one finds that 
when the following conditions are satisfied 27 

k= >> 1; k6M, << 1;1/k = 0(M)) (1) 


the equation for the perturbation velocity potential reduces to: 


‘ O*¢ ° 
M,? = V2 (2) 
of? 
where: 
WV = Mach Number at infinity ahead of the wing 

w = wing oscillation frequency 
h = half-chord of wing 
l = velocity infinitely far ahead of the wing (in posi- 


tive x-direction) 
k = reduced frequency 
6 = ratio of thickness or amplitude of wing to maxi 


oO = 
mum chord 
oe taal —_ 

Vo =(q — U)/U; q = gas velocity 

t, X, ¥, 3 = appropriately nondimensionalized time and space 


coordinates 


The boundary condition to be satisfied on the wing plan form is: 


o¢ t 
( *) = ikh(x, y)e* (3) 
Os J/-=0 


k(x, y) being the dimensionless amplitude of wing oscillation. 


The pressure coefficient takes on the form: 
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Cp = —2(d¢/22) 


It can be shown that the order of inagnitude of the terms in Eq 


(2) are as follows: 
O-¢ 
ot? 
] ¢ 
: 1 /1,*k? 
/;? Ox? - 
] 2 
1 /.17,74R2R? 
AI,? Oy? 
O*¢ 
] 
O22 


The derivatives with respect to x and y can therefore be neglected 
when the additional conditions: 
M,?k? >> 1 and V,27k? R2>> 1 6 


are satisfied. These requirements simply state that the wave 
length of the acoustic waves emanating from the wing must be 
small compared to all dimensions of the wing plan form. The 
equation then takes on the form 
we dy 7 
W2—* = 
ot- Os* 


with the boundary condition and the pressure coefficient given 
by Eq. (3) and Eq. (4), respectively 
An integral of Eq. (7) satisfying the boundary conditions can 


be written down immediately for ris 


The pressure coefficient is found to be: 


; 2Zikh(x, yye' Miz 
= Vu 
he corresponding expressions for arbitrary unsteady motion of 
the wing at points on the surface are 
g(x, y,0,t) = h(x, ye /M, 10 
Cp = 2ikh(x, ye /M, 1] 


The pressure difference between tip and bottom of the wing is 
twice this amount for the lifting case. Eq. (9) can be integrated 
to yield lift and moments 

Calculations have been carried out and checked against the 
two-dimensional theory of reference 6. Table 1 shows the re- 
sults at a typical supersonic Mach Number. The nomenclature 
is that of Garrick and Rubinow,® and subscript p refers to values 
obtained from the present simplified ‘‘piston theory.’’ As might 
be expected, ‘piston theory’’ proves more accurate the higher 
the Mach Number, but the order of accuracy shown in the table 


persists down to ., = 10/7 for k > 2. In view of the remark- 


able simplicity of these calculations compared with more exact 
theory (especially for oscillating finite wings in supersonic flow), 
the possibility of using them for applications like flutter analysis 


deserves careful consideration. Occasionally flutter speeds 


very sensitive to small charges in the phase angle of an 


are 
so that care must be 


aerodynamic quantity, like (14)’ + 7.1/,’) 
exercised in interpreting such applications 


REFERENCES 
1 Heaslet, M. A., and Lomax, H., Two-Dimensional Unsteady Lift Prob- 
Flight, NACA Report No. 945, 1949 
Two-Dimensional Rotational Flow at High Mach 
Journ. Mech. and Applied Math., 


lems in Supersont 
2 Goldsworthy, F. A 
Number Past Thin Airfoils, Quarterly 
Vol. 5, pp. 54-63, 1952 
Lighthill, M. J., Oscillating Airfoils at High Mach Number, Journal of the 
Vol. 20, No. 6, pp. 402-407, 1953. 


Aeronautical Sciences 
Quarterly Appl. Math., Vol. 5, 


4 Hayes, W. D., On Hypersonic Similitude 
pp. 105-106, 1947 





C 


tion 
prol 
(1 
Vth 
give 
bina 





READER 


n, C. C., Reissner, | rsien, H. S., On 7 Dime) ul N ad 
Motion of a Slender Body in a Compressible Flui Journ. Math. and Phys 
Vol. XXVII, No. 3, pp. 220-231, October, 1948 

Garrick, I. E., and Rubinow, S. [., Flutte ind O latin Lir-f 
t Foran A nl in Two-Dime tonal Super ic f NACA TN 
N 158, 1946 


In correspondence 
ir to Eqs. (1 4 


Prof. J. W. Miles points out he ha 


On the Nonsteady Climb of Turbojet Aircraft 


Angelo Miele 

Research Assistant Professor of Aeronautical Engineering 
Polytechnic Institute E 

June 10, 1954 





SYMBOLS 
speed of sound (fp 
drag coefticient at zero lift 
Le V wrreCpo 
I drag (Ib 


airplane elliciency factor 


ift coeflicient for maximum lift drag ratu 


Ne 


. Ta DO 


maximum lift drag rat 


acceleration of gravity s2.174 ft 
; height (ft 
altitude of the tropopaus« 35,332 ft 
he h+vV?'2 energy height (ft 
AM. Mach Number for best climb 
Ales = Mach Number for best climb determined from the equations 
of a steady flight 
l 2 
Mf, \ Mach Number for level flight with maximum 
1 ee 


CLep 
lift drag rati« 

I 1ir constant ™ 1717 ft eC R 

wing surtace (it 

thrust (Ib 

thrust at sea level (Ib 


thrust at the tropopause (Ib 


if velocity (fps 

Vv velocity for best climb (fps 

Ves velocity for best climb determined from tl equations of a 
steady flight (fp 

1 velocity corresponding to level flight with maximum lift/drag 


ratio (fps 


if weight of the aircraft (Ib 
i : 
ratio of minimum drag of level flight to thrust 
TI 
y ratio of speed for best climb determined from the equations of 


a nonsteady flight (1 to the speed for best climb determined 
from the equations of a steady flight (1 
derivative of the temperature of the air with respect to the 


iltitude 0.00385660°R ft for troposphere OPR ft for 
tratosphere 
ratio of specific heat at constant pressure to specific heat at 
constant volume 1.4 for air 
ispect ratio 
relative density of the air 
vy relative density of the air at the tropopaus 


f absolute density of the air (Ib sec? ft 


DISCUSSION 


ONSIDERABLI ittention has been devoted in recent years 
C to the study of the problems associated with the nonsta 
tionary climb of high-speed aircraft Some of the minimal 
problems investigated are the following: 

(1) The determination of the speed-height relationship |) = 
V(h) which minimizes the time necessary to fly an aircraft from a 
given combination of speed and altitude (1, / to another com 

Vo, h 

9 


(2) The determination of the speed-energy height relationship 


bination of speed and altitude ( 
= |(h,) which minimizes the time necessary to fly an aircraft 
from a given energy height level to another energy height 


level (h, with free choice of the terminal velocities, I, V 
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According to the terminology of the Calculus of Variations th 
first problem is a fixed end-points problem, while the second one, 
of lesser practical importance, is a free boundary value problem 
Mathematically analogous problems may be formulated for the 
climbing flight with minimum fuel consumption and for the steep 
est climb 

The writer! analyzed the fixed end-points problem by using 
transformation based on Green’s theorem rhis approach esset 
tially avoided the difficulties inherent to the use of direct or in 
direct variational methods. Both problems of absolute optim 
and conditioned optima were handled in this way 


K. J 


Salm 


In England 
Lush,* following a previous study by Kaiser,® treated thi 


fixed end-points problem with a graphic-analytic metho 
based on the concept of energy height (A, He also analyze 
the free boundary value problem. Later the writer* ® general 
ized his previous solution to the climbing flight of minimum fuel 
consumption, to the stecpest clinb and to some speci il type ol 
Additional contributions to the san 


Kelly® and more r¢ 


nonsteady gliding flight 1 


of problems were developed by 
and Garbell.* 


The following hypotheses are used: the small angle betwee: 


category 


cently by Rutowski 


the vectors thrust and velocity is neglected; the weight of the 
aircraft is assumed constant; the curvatures and the squares of 
the path inclination are considered negligible with regard to their 
effect on that part of the drag depending on the lift; the thrust 
is assumed a function 7 7(V, h) of the velocity and altitude 
only the aerodynamic lag is disregarded; only flight paths 
contained in a vertical plane are considered 

The distribution of velocities V V(h) which, under the afore 
mentioned assumptions, minimizes the time integral in the fixed 
end-points problem is composed of three branches Pwo of 
them, the initial and the final branches, depend on the boundary 
conditions of the problem The third one, the central branch 


is independent of the b.c., being defined! by 


re) Gan _Vfo (7 . ory | 
OV iW , g | Oh fl | 


With a transformation of coordinates from the (VV, h 
the ( V, Ht, 


plane to 
plane, Eq. (1 


[oO [= ~ PY) 
() » 
Lor i J r 


For the free boundary value problem the distribution of v« 


may be shown equivalent? to 


locities V V(A,) minimizing the time integral is composed by 


a unique branch defined by Eqs. (1) or (2 


Due to the lack of simple analytical expressions for the thrust 


and for the dr ig, a gener il solution of Eqs 1) or (2) has not 
yet been found. Graphical methods* * are to be employed for 


each particular aircraft, above all when compressibility effects 


ire of importance, as is the case for turbo-jet aircraft at the 
altitudes in the neighborhood of the ceiling 

On the contrary a rather simple solution, of use for preliminary 
design calculations, may be obtained if the following engineering 
issumptions are accepted for the thrust and the drag 
iltituce 


(1) The thrust is assumed to depend only on the 


cording to the following relationship 


T = Ko’ (3 
where 
K = 7), r = 0.75 for tropospheric flight 
r Ts 
K = , » = 1 for stratospheric flight 
CO x 
» 


(2) The polar of the aircraft is supposed parabolic with con 


stant coefficients Cyp, « This assumption leads to the following 
expression for the aerodynamic drag: 
I 2W? 


D = , Coo pSV T 9 (4 
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Fic. 1. Mach Number (.\/,.) for best climb (troposphere) Use the lower Jf, scale to determine /. = f(x, \/ Use the 
upper ./, scale to determine } = f(x, AV/,). 


From Eggs. (1), (3), (4) and from the Equations of the standard 
atmosphere, the speed for best climb (V,) can now be calculated. 
Results may be expressed in the following simple nondimensional 
form by introducing the Mach Number (./,) in place of the ve- 
locity as a variable and by using the parameters x, .\/, already 
defined in the list of symbols: 


M6 +AM,.4 + BM,.? + C=0 (5) 
where: 
3 M,? 
1 = — 2r 
ym x 
2 M,? 
B= — JM,‘ 
ym Xx 
(6) 
M.‘ 
C=- 
ym 
aR 
m = 


Eq. (5) is a cubic equation in .1/,?, which may be analytically 


solved with Cardano’s formula or with equivalent trigonometric 


expressions. The coefficients B and C are always negative 
For A(discriminant) < 0 Eq. (5) has two real negative roots and 
For A 


jugate roots and one real positive root. 


one real positive root. > 0 Eq. (5) has two complex con 
The only positive real 
root of Eq. (5) has been calculated and plotted in graphs 1 and 2 
which enable one to carry out a rapid calculation of the optimum 
climbing program 

The Mach Number for best climb (.1/,) depends on two param 
the ratio of minimum drag to thrust (x) and the 


with maximum lift/drag ratio 


eters only: 
Mach Number for level flight 
(.\7,). These two parameters are fundamental in the study of 
the mechanics of flight of jet-propelled aircraft, since they em- 
body the effects associated with the shape of the aircraft (Cpo, 
\, e), the wing loading (W/S), the specific thrust (7 1”), and the 
properties of the atmosphere in which the aircraft is flying (a, p 
It is of interest to compare the actual results with those pre 
viously derived from the old quasi-uniform analysis of the climb- 
for best 


ing flight. In the quasi-steady approach the speed ( V’;. 


climb is defined by: 
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Fic. 2. Mach Number (.\/,) for best climb (stratosphere.) 

Use the lower \/, scale to determine 1/7. = f(x, VW, Use the 
upper .l/, scale to determine }° = f(x, M 

fe) (’ Vv- =~) = 

=s ¢) (4 
oV W h 

For a parabolic polar and for (O7/OV), = 0 the velocity (V5) or 


the Mach Number (.\/,,) solving Eq. (7) is given by 


Ves Mes i+ ¥ 


ie Pek 3x =) 


fe, being a M, and AM, being a function of x, 
M, it follows that the ratio Y = (M,./Mes) = (Ve/ Ves) 
Y= f(x, \/,) has been plotted 


function of x, 
isalsoa 
function of x, \/,. The function 
in the graphs 1 and 2, which show that the speed for best climb 
determined with the actual nonsteady analysis is higher than the 
one determined with the old quasi-steady treatment. The rela- 
increases with the altitude up to 


1) and 16-20% for 


tive difference (V. — Ves)/Ves 
about 8-10% for tropospheric flight (Fig 
stratospheric flight (Fig. 2) 

For a practical use of the present results the following remarks 
are of interest: 

(1) Ath = hx (tropopause) there are two velocities solving Eq 
(1 One (Ve trop) 
flight, the other one (V¢ strat) is the initial velocity of the strato 
In Ref 
spheric climb and preceding the stratospheric climb, the aircraft 


is the terminal velocity of the tropospheric 


spheric flight 4 it was shown that following the tropo 


is to be accelerated at h = he from V¢ trop to Ve strat in order to 
achieve the optimum climb performance 

(2) The solution of Eq. (1) is much more sensitive to the /oca/ 
value of the thrust (7) than to its derivative (d7/dh) with respect 
tothe altitude. Therefore, even if for a given turbo-jet engine the 


exponent r of Eq. (3) is different from the assumed value r = 


and (2) may still be used with good approxima 


(W/TEmar.) be determined 


0.75, graphs (1 
tion, provided the parameter x = 
on the basis of the effective thrust 7 at altitude from the graphs 
supplied by the turbo-jet maker 

(3) Results 
speeds, because the hypotheses OCpo/OM = 0, 0e/0\J = 0 are 


are to be handled with caution at high subsonic 


no longer satisfied. When compressibility effects are of engi- 
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neering importance graphical methods of solution of Eqs. (1) or 
(2) are to be preferred to the present one 
4) For small values of (O7'/OV), the 
with the velocity may be approximately accounted for by de 
= W/TEma:z 


Eq. 5 or graphs 1 and 2) with trial and error method 
1 gra] 


variability of the thrust 


termining the parameter x ind the optimum ve 


locity V, 
When, on the contrary, the thrust is a rapidly varying function 
of the velocity, graphical methods of solution of Eqs. (1) or (2 
are to be preferred to the present approach which neglects the 
terms containing (O07 /OV ), in the equation defining the optimum 


speed 
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perturbation velocity po 
ta 


HE INTEGRAL EQUATION for the 


tential ¢ which satisfies the condition of tangential flow 


point Q on the surface of a wing is given by Evvard.' 


* A complete discussion of this method will be published soon in the ‘*Pro 
ceedings of the Second Canadian Symposium on Aerodynamics,"’ published 


by the Institute of Aerophysics, University of Toronto 
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The region of integration is indicated in Fig. 1 If the forward 
running Mach lines from Q are reflected off the wing leading 
edges, as shown in Fig. 2, the wing surface is divided into a 
series of oblique rectangular regions, which are denoted by R,, 
Ro, R3, ete. 


of these regions is calculated by means of Eq. (1 


The contribution to the potential at Q due to three 
and a knowledge 
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Fic. 3. Comparison of approximate and exact solutions for a 


flat delta wing. 


of the strength of the upwash field between the Mach cone and the 


leading edge of the wing 


tain 

la dr ¢ 
irs ) 
ade = ~ 

Vor V on = 

dg. = 0 > 
_ ] 4 
ag — oe + 


rhus, the integral equation for ¢ at Q may be written 


. . 
+ | de : { d ¢. T 5 
JR JR 


Phe magnitude of the terms in this series decreases rapidly 
due to the presence of the radical in the denominators, and to 
the decreasing size of the more remote regions 

rhe approximation made by Evvard? in connection with vawed 


wings Was to retain only the first term of the series. This gives 








Fic. 5. 


Subdivision of region 3. 





see Eqs. 38 and 39 of reference 1) to ob- 














READERS’ 


of the but is con- 


good results over useful range parameters, 
siderably in error for highly swept wings, or for wings at low 
upersonic speeds rhe second approximation, in which the 


ul 
integration is carried out over both Regions 1 and 3 is shown to 
ive a significant improvement to the results for flat delta wings 
t small values of the parameter 8 tan 8, as is illustrated 
rhe 


faces by 


in Fig. 3 


method may be extended to wings with nonplanar sur 


further subdividing the regions as indicated in Figs 


tand 5. For example, Region 1 is subdivided into WV equal 


pal illelograms of sides a, 6; the size being such that the local 
ngle of attack a , may be considered constant in each Inte- 
grating Eq. (2) over one of these small parallelograms 
$Uam, 
Ago = \ ey ae ee TY aS Yes 
eM 
6) 
Noting that, for wings with straight subsonic leading edges 
V : 
= v — Bky 
3VU(1 +k 
V/ 
= GRY 
3M R 
whcre 
l » 
k= =V iV, -—1 


8 tan 6 


y) are the coordinates of Q, 


¢ = : VV X* — (BR )*¥* Sm, Ss 
: j7r(1 R 
where 
VV - \ kw 
Ee Ae ae ae Lae = 28 ‘ 
I l V NVM 
rhe pressure coefficient at Q due to Region 1 may now be cal 
culated 
2 ¢ s 
ae = — . 
l Ox 3r(1 +R 
S oe 7 
Vl —(k x 10 
V1 — (hk Ox 
here = +V \ 


Similarly, the pressure coetiicient at Q due to Region 3 is given 


VI — (ep Ox 
Vhere 
1 P Oo 
= = = 
ro = — 


1+kt Kp 1 — kt 
V 7 1 — kt Py (a ’ 1 + ki 


The method may equally well be applied to nonplanar wings 


with curved subsonic leading edges by laying out the regions on a 
planform drawing of the wing and measuring off the required 
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quantities directly rhe pressure coefficient is obtained by a 


final numerical differentiation of the potential so calculated 


In a practical computation of the pressure distribution over a 


cambered wing, a drawing of the wing planform is divided into 


Regions 1 and 3, and these regions further subdivided to give 


the required net If the local angle of attack contours are plotted 


on the wing drawing, the angle of attack matrices a», ,, @», ¢ May 


be read off directiy The summations S and ire then 


computed, differentiated numerically, and placed in formulas 


8) and (11), respectively, to give the desired pressure coefficients 


It has been found that for most problems very little error is intro- 





duced by using a coarse net, although ./ or N should not be less 
than 3 
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Probe Interference in Measurements in 
Supersonic Laminar Boundary Layers* 


M. V. Morkovin and W. S. Bradfield 
oO 


Department of Aeronautics, The Johns Hopkins University, and 
Rosemount Research Center he University of Minnesota 
Respectively 

June 15, 1954 


Rees, EFFECTS OF THE SIZE of the total-pressure probe 
used for measuring velocity and other profiles of super 
ceived considerable attention 


the 


sonic laminar boundary lavers r¢ 


e.g., in references 1, 2, and 3 Inference of existence and 


partial evaluation of these effects were based (a) on comparison 


of results from different size probes in references 1 and 2 and (b 
on compari:on of results with boundary-layer theory (which 
cannot count for the complicating leading-edge effects) in 
references 1, 2, and 3 rhe purpose of this note is to present 
direct visual evidence of strong interference near the wall and dis 
cuss its implications 

Briefly, local, three-dimensional separation under various 


as observed in the Johns 


Ms= 1.72 


probes (hot-wire and total pressure) w 
> by the 


Hopkins University tunnel over a flat plate at 


light-oil lamp-black technique These observations were re 
peated and carefully documented for even the smallest practical 
total pressure probes’ at the University of Minnesota at 1/ 

Ww illy supported by the Office of > R \ 

the Office of Scientific Research, USAI 


— — 





Fic. 1 Pattern 


for probe tip t mils above flat plate 
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Fic. 3. Probe without background pattern. 
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Fic. 4. Velocity profile inferred from pitot readings and oil film 
behavior at positions 5, 9, and 16 mils above cone surface 
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3.05 with the aid of a short focus schlieren system and the chiig 


film technique. Fig. 1 gives the top view of the probe and thy 
pattern when the center of the probe, 0.008 in. high and 0.020 in 
was 0.004 in 


wide (see photomicrograph in reference 3), above 


the plate as determined by a 20X microscope during running 
In Fig. 2, the distance from wall is 0.007 in. and Fig. 3 shows the 
probe without the background pattern Visual observations 
during running, accompanied by vertical and horizontal motion 
of the probes identify the ‘‘island’’ pattern as one of local three 
dimensional separation with back flow immediately under the 

The boundary 
and the Reynolds 


Number based on distance from the leading edge was 0.97 


probe, extending upstream of the probe tip 


layer thickness in these cases was 0.025 in 


10°. Fig. 4 displays the velocity profile of a thinner boundary 


layer measured with the same probe on a cone at MW = 3.05 and 
diagrams illustrating the behavior of a thin oil film on the under 
side of the pitot tube. As the probe approached the surface the 
oil film initially distributed as shown in the top sketch began to 
collect, toward the leading edge of the under surface of the tube 


as shown, indicating reverse flow at the probe as well 


In case of larger probes upstream separation tends to persist 
even for positions of the probe nearer the edge of the boundary 
laver. From similar experiences with hot-wire probes one can 
surmise that the velocity profiles from the thinner layers near a 
leading edge reported in reference 3 as ‘‘transitional” (p. 8) and in 
reference 4 were actually obtained in a /aminar layer with the 
extent of separation varying as the vertical position of the probe 
was changed. Even such more extensive separation is difficult 
to observe with schlieren systems because its lateral extent is 
small in comparison with the width of the plate. This is prob 
ably the reason why this effect has not been reported heretofore 


It has thus been demonstrated that when a probe of given 
height / is located at a distance from the wall of the order of h, 
a strong interference including local separation and back flow 
This finding supports the disregard of pitot tube data 


Pitot tube data 


occurs 
near the wall as practiced in references 1 and 2 
may still be subject to error when such separation does not occur 
or occurs downstream of the probe tip. Lighthill’s discussion of 
viscous two-dimensional upstream influence without separation’ 
can perhaps be used as a guide keeping in mind, (a) that three 
dimensional effects should be milder, (b) that the pressure dis 
turbances in front of the probes correspond to large x and y gra 
dients to which boundary layer theory is not strictly applicable 
Lighthill’s upstream influence dies down exponentially, extend- 
ing a distance comparable to ten boundary-layer thicknesses, 
depending upon Reynolds Number and Mach Number. Even 
if the effect were considerably decreased by the three-dimensional 
relieving action, some thickening of the layer with consequent 
oblique pressure waves in the supersonic part of the boundary 
layer should result. Thus, apart from the corresponding dis 
tortion of the layer the assumption of constant static pressure 
throughout, used in transforming the pitot data into velocity 
and mass flux information, is affected 


Fig. 4 of reference 2 and Fig. 22 of reference 1 show an ap- 


parent ‘‘overshoot’’ of momentum flux and velocity above the 
free stream values for the relatively larger (0.008 in.) probes 
even when completely outside of the boundary layer (6 ~ 0.04 in 

(A similar observation was made in reference 4.) Since (a) this 
effect is not noticeable for the smallest probes in the same bound- 
ary layer and since (b) it vanishes as the larger probe moves still 
further away from the wall, the effect must be ascribed to up- 
stream influence through the boundary layer. The pressure 
waves previously mentioned and the consequent small decrease 
in local free stream Mach Number would account for the higher 
pitot reading and the fictitious velocity overshoot. Using the 
Rayleigh formula, one can show that the error in indicated 
Mach Number, A.V, is related to the rise in static pressure, Ap;, in 


front of the pitot mouth by the formula: 
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, similar formula in terms of the pitot tube reading shows that the 
overshoots indicated in Fig. 4 of reference 2 could be explained by 
1-5 


: per cent rise in static pressure p, in front of the pitot 
mouth 


The preceding conjecture as well as the general up 


stream effect can be checked by taking simultaneous readings of 
pitot pressure at each y position and of static pressure from wall 
wifices under and slightly upstream of the pitot tip 

If there is a demonstrable effect of upstream influence when a 
relatively large probe is outside of the layer, it would follow that 
similar or even larger effects could occur throughout the boundary 
laver and that some effect should be present for smaller probes 
The presence of this effect and that due to the leading edge makes 
the near agreement with theory (which takes account of neither 
in references 1, 2, and 3 the more remarkable 

rhe question remains as to the non-dimensional parameters 
properly characterizing the above interference, with or without 
separation. Clearly, the probe height / has to be non-dimension 
ilized with respect to some characteristic layer thickness in order 
to account for the observed effect of ‘‘relative size Just how 
this should be accomplished is not clear since it involves ques 
tions of pressure disturbance strength as function of y and the 
upstream viscous influence. The width of the probe must be 
significant since its increase means closer approach to the ob 
two-dimensional critical disturbance for sepa 


served pressure 


ration, 


Ap =4op V 


(or equivalent), which is readily transgressed. Cases have been 
observed where tendency to local separation near the probe was 
reinforced by pressure disturbances from the bulkier supports of 
the probe as much as 30 boundary layers downstream from the 
tip. Such spurious effects would, of course, spoil correlation of 
data based on tip size alone 

he framework of Lighthill's work indicates dependence on 
Reynolds Number R,, Mach Number .V/, and wall-to-free-stream 
lr; 


s, through which a given (linearized) pressure disturbance at wall 


temperature ratio 7), for instance, the upstream distance, 


diminishes by a factor e~!, measured with boundary layer thick 
ness as a unit, depends on the parameters approximately as fol- 


lows: 
1.3 oo ; M+ 2 To \° 
a/ Mi — 1 11 20M -1\ % 


There does not appear to be data over a sufficient range of Mach 
and Reynolds Numbers to validate the simpler correlations at 
tempted so far 

In the case of turbulent boundary layers the relative size of 
probes to boundary layer thickness is usually small Further 
more, the upstream influence, in units of boundary layer thick 
Still, the effect 


probes in super- 


ness, is also known to be very much smaller 
of the sharp pressure gradients generated by 
sonic layers is likely to alter the simple field postulated in evalu 


ation of the readings. Velocity profiles obtained from probes of 





»bes 
in.) 
this 
nd- 
till 
up- 
ure 


her 
the 
ted 


, in 





different sizes in reference 8 indicate that there is a significant 
effect of size for probe-diameter to layer-thickness ratio larger 
than 1/20 at M = 2 
0.6 were available, the range of discrepancies might be found more 


If data for velocity ratios (u/V,, ) less than 


severe 
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Electrical Resistance Networks for Beam and 
Column Problems 


R. H. Scanlan 

Centre National de la Recherche 
Bagneux, Seine, France 

June 21, 1954 


cientifique, Chatillon-sou 


M** ARD AND BOSCHER, in a series of condensed notes, 
have exhibited the essential characteristics of a useful all 
resistance electrical network which represents various two 


dimensional fields. The fields may represent the potentials for 
fluid flow or, 


Their basic idea will first be briefly described 


stress or for example, the deflection pattern of a 
plate or membrane 
here 

Consider two two-dimensional all-resistance networks I and II 
superposed and tied together at their nodes by other resistances 
(Fig. 1). Let us permit the entrance of outside currents from a 
source at potential U’ through still other resistances (® and R 
above and below the two networks. One node of the system ts 
shown in the figure 

If ®, R, R are very much larger than r and p, it is possible (and 
this is currently being exploited in several directions) to solve, 


using the form of finite differences, two simultaneous systems 


governed by the equations 


o ol re) Ol - 
m + n = AV + F 
Ox OX Oy Oy 
re) ol 0 Ol 1 : 
p 4 Q = —-(C) + BV +G 
Ox Ox oy oy 


where m,n, p, g, A, B, C, F, Gare known functions of x and y and 


V, V are the node voltages of the lower and upper networks, re 
spectively, which correspond to the values of the functions 
sought From a general point of view, and one which aids 


physically in visualizing the nature of these networks, it can bc 
network II 
regulates the nature of the voltage 
voltage network”’ I, it being assumed that the lower network has 
Typical values of the voltages 


said that is the “higher voltage network’’ which 


distribution in the “lower 


practically no effect on the upper 


and resistances used make this clear: r = 1,000 ohms, R = 1 
megohm, p = 100 ohms, ® = 100,000 ohms, R 1 megohm 
Voltages used may be U’ = 100, F about 3, V about 0.01 on the 


scale of maximum voltage taken as 100 
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Here we note in particular the use of a single line of nodes cut 


out of one of these two-dimensional networks. This system (sec 


Fig. 2) solves the classical beam equations (where /: is deflection, 


M moment, p load, ete 


d?\ 

dx? i a 
d*h 

El = M(x 
dx? 


Fig 2 shows a setup used to obtain the influence coeflicients of an 
aircraft wing spar. The accuracy obtained is of the order of 1 per 


cent. If the voltage l’ is called 100 (actually, as used, in the 
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neighborhood of 30 volts), the voltage representing the beam d 


flection # can be read to within 0.001. Low frequency (120 evel 
a.c. electrical feed permitting the use of a cathode-ray oscillograp] 


as a zero instrument is very convenient for this type of network 


the low frequency reducing the possibility of a.c. effects 
essentially d.c. network 
It is easy, with a sufficiently large number of network nodes, t 


obtain the vibration modes of a beam by an iteration schen 


which resembles ver\ closely the classical matrix iteratioy 


Through proper choice of the resistances R it is also possible t 


bring in easily the effects of rotary beam inertia. Naturally, vari 


able section rigidity EJ causes no more difficulty than constant 


rigidity, since R is chosen proportional to EJ at each node. Thy 
typical beam boundary conditions are regulated simply by po 
tentiometers. Higher modes are amenable to ‘‘sweeping matri 
type of operations 

Another resistance network of a simple nature is useful in repre 


senting the torsion equations of a bar. This is shown in Fig. 3 for 


a cantilever bar having four discrete sections of different stiffnesses 


GJ 
The voltage | corresponds to the angle of torsion 6 of the bar 
the currents which enter (=l°/R 
torque 7’, and ¢ corresponds to /,G/. Here also we take R >? 
This network is good for torsional vibration problems, worked out 

by iteration as before. 
If the idea of an elastic or torsion axis can be tolerated in a giver 
case, simultaneous iteration with the two networks like those of 
) 


Fig. 2 and Fig. 3 and the standard equations of coupled vibration 


d*h 
d? avs] BI ( F ‘)] = w(ph + SO 
ax a 


GJ(d6/dx) = —w Sh + Je 
an 


will yield rapidly the coupled mode of the system. (yp, S, and J are, 
respectively, the mass, static moment, and moment of inertia per 
unit length 

Other characteristic value problems of structures can be at- 
tacked with the beam flexural network.* Among these, also solved 
by rapid iteration, are: critical loads of columns with any rigidi- 
ties, end conditions, and distribution of longitudinal loading (the 
resting upon elastic 


bar may also be assumed vibrating or 


lateral supports of any stiffness); the buckling of arches in their 


planes; the lateral buckling of thin beams. Cases of beam-column 
deflection can be accurately solved by iterative approximation 
In view of the somewhat complicated form occasionally taken 
by electrical networks it is of interest to note that certain results 
lesser 


of fairly wide application can be obtained, as here, by 


means 
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Rin Megohms; 7 and p in ohms; and ® = 10,000 ohms. 
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Errata—‘‘ Differentiation of Experimental Data 
by Means of Higher Order Finite-Difference 
Formulas’’ 





Chi-Teh Wang and Daniel F. De Santo 

Daniel Guggenheim School of Aeronautics, New York University, 
New York 

June 23, 1954 


N THEIR READERS’ FORUM NOTE,' Wang and De Santo pre 
sented a method for differentiating experimental data 
using the method of higher order finite-difference approximations 


In that note they employed as an example the function 


by 


x = 6 sin? + 0.5 sin 3t (0 is 3 l 


It was subsequently pointed out by Schurmann and Smith? that 
the values of dx/dt and d?x/dt? calculated by the method of refer 
ence 1 at the first, second, next to the last, and last points on the 
respective curves might not be correct. The authors 
have since recomputed the derivatives at these points and found 


present 


that numerical errors had indeed been made in calculating these 


values. The correct values are presented in Fig. 1 It may be 
seen that the accuracy is now good over the entire range, includ 
ing the points in the vicinity of the ends, and that the resuits 


compare favorably with those given by the method of reference 2 
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A Note on Rolling Pull-Out Calculations 


George S. Campbell 
California Institute of Technology 
June 26, 1954 


| 


Pasadena, Calif 


PAPER ENTITLED ‘‘Estimationu of the Maximum Angle of 


NA 
I Sideslip for Determination of Vertical-Tail Loads in Rolling 


Maneuvers,”'! it was shown that the three linearized lateral 


equations of motion, including product of inertia terms, may 
not provide a sufficiently close estimate of the sideslip angles « 

pected in rolling pull-outs rhe results presented in the present 
note indicate that the nonlinear, four-degree of freedom system 
of 


possibly just 


a standard comparison in reference | 


of equations used as 


are as questionable as 


may provide results that 
those obtained from the linearized equations 

The exact equations of motion representing the six-degrees of 
freedom of a rigid aircraft are conveniently given in a recent text 
by W. J The yawing-moment equation, for example, 
has the following form under the assumption of body axes which 


Duncan.? 


are principal axes of inertia 


V =/;i i & 


) pg (] 


When the complete expressions for the angular velocities p, 


and r, as given in reference 2, are introduced into Eq. (1), the 


yawing-moment equation becomes 
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Effect of nonlinear terms on the aerodynamic yawing 
moments during a rolling pull-out 


Fic. 1. 


N = IAW cos ¢ cos 0 6 sin @ — po sin @ cos 8) + 
(lr — 1,)W sin @ cos AW sin @ o) 


(iz + I, 1,)o6 COS @ (I, I, + 1,) Ve cos @sin@ (2 


For comparison, the yaw-equation used for the nonlinear 


calculations of reference | is 


N = IW — (I; — 1,)66 (: 


This form is actually subject to some serious linearizing approxi 


nations, for it is obtained from Eq. (1) by assuming 


p =| 
gq = 6 
r = pl 
Equations (4) are valid approximations only if ¢, 4, and Ware small 
quantities, and it is evident that ¢, for example, may reach large 
values in rolling maneuvers 

Finally, the usual linearized equation of yawing motion may 


be written 


In order to determine whether the nonlinear equat’ons used in 


reference | are sufficiently close approximations to the complete 


2 _ 
Pyy = 


Summing the drag forces and equating to zero, 


ay: Scie 


' (Cr) 
aSCp 
dt 
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six-degree-of-freedom equations of motion, the results for thy 
case of a rolling pull-out presented in Fig. 4 of reference | may 
1 of the present note repre 


be used. Curves presented in Fig 


sent the aerodynamic yawing moments required to produce th 


specified rolling pull-out maneuver These yawing moments 
are calculated from the three expressions |Eqs. (2), (3), and (5 
for the inertial reaction during the maneuver. It is evident 


from these results that the terms which are neglected when Eq 
(3) is used for time history calculations are not necessarily negli 
gible for the case of a rolling pull-out For if these terms actually 
had no significant effect on the time histories, the values of th 


yvawing moment as calculated from Eq. (2) 


would coincide with 
the curve obtained from Eq. (3) 

It is further seen from Fig. 1 that the use of only one nonlinear 
term, as in reference 1, may produce results which could be in 
greater disagreement with results obtained from the complet 
nonlinear expression than would have been obtained from the 
much simpler linear equations 

From these results, it is suggested that if one is willing to accept 
the complications arising from the inclusion of any nonlinear 
terms in the equations of motion, then the use of a// the nonlinear 
terms is probably necessary in order to obtain reliable results 
The question of whether four degrees of freedom are adequate 
for the rolling pull-out is a separate consideration beyond the 


scope of this note 
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Airplane Stopping Distance 


B. Saelman 
Engineering Consultant, North Hollywood, Calif. 
June 29, 1954 


_—— DRAG FORCES acting on the airplane during the decelera 
tion from touchdown to the stop are: 

(a) Rolling resistance of the nose gear tire 

(b) Rolling and brake resistance of the main wheels (skidding 
force if brakes are locked) 

(c) Aerodynamic drag on the airplane 


d) Positive or negative propeller or jet thrust 
In reference 1 a formula for stopping distance is derived base« 
on the assumption that the braking force is given by 


F = pW — L) 


This will apply to airplanes with brakes on both the main and 
nose wheels. However, when the brakes are applied only to the 
main wheels, as is usually the case, this assumption is not real 
istic and, if used, can lead to considerable error. For this case 
formulas are derived below, and have yielded answers which have 
closely checked measured stopping distances when an average 
cocthcient of friction of 0.3 was used 

Taking moments about the nose wheel and solving for the main 


wheel reaction, 
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n Eq In these equations the terms are as defined in reference 1 and in Fig. 1. Letting y = dx/dt, Eq. (2) can be simplified to 
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Assuming 8 positive, integrating Eq. (3), and using the initial con- If p 0 
dition that dx/dt = V; whent = 0, | 
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ies dt H G tanh of : — 
a where 6; = 8 is substituted in Eq. (8) 
sing > termi . iti = () ni-=_ . ‘ , , , . . 
[ ing the terminal condition that dx /dt when ! 1, we have Integrating Eq. (6) and using the terminal conditions the 
the time of the ground run stopping distance is: 
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In reference | and in this discussion it was assumed that yu Fig. 2 is taken. Fig. 2 shows the variation of torque with speed 
vas constant and independent of velocity, and that ~ was con in a typical dynamometer test wherein brake pressure is held at 
istent with the ‘‘torqueability.”’ of the brakes. This latter phe i constant value throughout the stop rhe solid line shows 
nomenon is quite significant and is often overlooked. It is well the variation in torque (at constant ike pressure) which 


deseribed in an unpublished work by K. S. Carter from which — results from the combined effects of brake dise speed chang 
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The 
dotted line shows the torque which is usable on the airplane, 
and is defined as the maximum torque which can be applied to 
the brakes, at each instant of the stopping period, without 
This torque is a function of airplane 


and temperature changes which occur during the stop 


danger of skidding the tires 
speed and wing lift at the instant, together with the coefficient of 
friction which can be developed by the tire against the runway 


surface 


REFERENCE 


Reverse Pitch Propellers as 
Vol. 12 


Crockett, Harold B., and Bonney, Arthur E 
Airplane Landing Brakes, 
No. 4, p. 441, October, 1945 


Journal of the Aeronautical Sciences 


A Reply—Load Distributions on Delta Wings 
at Supersonic Speeds 


Julian H. Kainer 

Aerodynamics Engineer, Ryan Aeronautical Company, San Diego, 
Calif. 

July 2, 1954 

ie A RECENT PAPER! in the Readers’ Forum, Mr. Zienkiewicz 
derived the load distributions on some twisted and cambered 

delta wings with supersonic leading and trailing edges. Upon 

comparing the results of Eq. (8) therein, Mr. Zienkiewicz errone 

ously concluded that his results differed from the load distribu 


tions we presented previously? by a factor of 1/8. Both results’ ? 
were compared to the load distribution for the classical Ackeret 
region (Ap/qg = 4a/8) from which comparison it was shown that 


our results? reduced to Ap/q = 4a 


AERONAUTICAL 


SCIENCES NOVEMBER, 

A similar observation was first brought to my attention by Dr, 
Godfrey Lance who has completed some comprehensive studies 
which include the effects of both symmetric and antisymmetric 
downwash distributions for twisted and cambered wings having 
Through 
it was established that our results 


subsonic, transonic, and supersonic leading edges 
correspondence with Dr. Lance, 
are correct 

It is the purpose now to emphasize the fact that we defined for 
loading function and abbreviated it 


The ‘effective’ function (i.e, 


convenience an “‘effective”’ 
with an asterisk as a superscript 
pressure distribution, velocity potential, horizontal perturbation 


velocity, section properties, and overall wing characteristics ) is de- 


: ; ‘ : ; Ap\* 
fined as the product of 8 and the desired function, i.e 


g 
1 

kiewicz should have been made was that (Ap/q)* = 4a@ which is 
the same as Ap/g = 4a/8! To indicate this more concretely let us 


look at a few specific formulas from our results:* the lifting pres- 


Ap —" _ : 
8 (See Symbols? Thus the observation Mr. Zien- 


sure distribution was given as 
12a) 


and the generalized horizontal perturbation velocity parameter 


T , 
V" (B/x (dy a) 


For the classical limiting case mentioned above, the conical func- 


was defined as 


Bit) (13b) 


tions, B,(t) were indicated by the barred quantities, B;(t). Thus, 


for 7 = 0, Eq. (13b) becontes 


Substitution of which into Eq. (12a) gives the expected result, 
Ap/q = 4a/8. This result was graphically demonstrated? in Fig. 


» 


In the original version,* the relationship between the asterisked 
function and the corresponding function without asterisk was 
more explicitly defined than in the more current version.? In the 
Symbols of the former, the asterisked loading function was de- 
scribed as an effective loading function defined as the product of 8 
and the loading function in question. It is hoped that the omis- 
sion of this definition will not continue to be a source of con- 


fusion. 
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